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Introduction 



Many-body physics is a branch of physics whose scope is to understand 
physical phenomena where a number of interacting bodies is present. The 
presence of the interaction is what makes the description of such systems 
challenging but at the same time exciting. Interacting particles can give 
birth to new physical process which cannot be simply described as the sum 
of the behaviour of each single element. The superconducting phase at low 
temperatures, plasmon peaks in the absorption spectrum, Mott transitions 
are only a few examples. New physics emerge as a result of the coherent 
behaviour of the many-body system. 

The description of interacting particles requires sophisticated many-body 
techniques and the exact mathematical solution to the problem is almost 
never available: approximations are needed. To construct a practical ap- 
proximation one need to have some clue as to which are the most relevant 
phenomena, which are the physical aspects that can be discarded and which 
can be treated in an approximate way as perturbations. Often even the ap- 
proximate equation cannot be solved analytically, a computational approach 
is needed. Computational Physics can be seen as an approach which stands 
in the middle between Theoretical Physics and Experimental Physics. Some 
of the results presented in the present work have been obtained through a 
numerical approach. 

In the present work we will describe some of these techniques with a 
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focus on a specific phenomenon: tlie description of double excitations in the 
absorption spectrum. Double excitations are a peculiar effect of interacting 
systems which does not have a counterpart in non-interacting ones. The 
optical absorption spectrum of a system is obtained by shining light on it. 
At the microscopic level photons hit the electrons which sit in the ground 
state and change their configuration. If the light source is not too intense 
this can be described in linear response; that is only "one photon" processes 
are involved, only one electron per time can be influenced. Here is where the 
interaction comes in. The hit electron is linked to the others and so other 
process take place, one of these is the appearance of multiple excitations. 
These are, strictly speaking, virtual processes as the real time evolution of the 
system is different from the one described. Nevertheless the physical effect 
is there and can be measured as extra peaks in the absorption spectrum. 

Double excitations is not the only subject of the present work and under- 
standing is not the only scope of many-body physics. The same techniques 
can be used to make accurate quantitative predictions of the behaviour of 
a material. These allow us to control physical phenomena and possibly to 
use them in technological applications. In the second part of the Thesis we 
focus on the application of more standard techniques to the description of 
carbon nanotubes (CNTs). In particular we focus on the effects of magnetic 
fields on CNTs. 

CNTs are quasi ID-systems composed by carbon atoms which have been 
discovered in They have the shape of a hollow cylinder with a nano- 

metric diameter (10~^ m), a micrometric length (10~^ m) and the thickness 
of a single atomic layej^. What makes such objects so interesting is that they 

large percentage of academic and popular literature attributes their discovery to 
Sumio lijima of NEC in 1991 [HE], however already in 1952 L. V. Radushkevich and V. 
M. Lukyanovich published clear images of 50 nanometer diameter tubes made of carbon 
in the Soviet Journal of Physical Chemistry, the publication however was in Russian. For 
a detailed review on the discovery of CNTs we address the reader to Ref. [3]. 

^Here we refer to single-walled CNTs. Multi-walled CNTs, which are composed by 
concentric singlc-wallcd CNTs exist too. 
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are mechanically very strong and stable. These properties makes them ideal 
system both for many possible technological approaches and for testing the 
physical behavior of electrons in ID system as well as in cylindrical topolo- 
gies. In this work we arc in particular interested in the effect of magnetic 
fields related to topology. 

Under the effect of a magnetic field electrons delocalized on a cylindrical 
surface display a peculiar behaviour, known as Aharonov-Bohm effect. The 
Aharonov-Bohm is a pure quantum mechanical effect which does not have 
any counterpart in classical physics. In CNTs the Aharonov-Bohm modify 
the electronic gap and so can be used to tune the electronic properties. 
Though a model able to account for such process is available in the literature, 
in the present work we will describe the effect of magnetic fields "ab initio " . 
Ab initio is any approach which describes the physics starting from first 
principles and without the use of any external parameter. As pointed out 
in the first part of the introduction the exact solution to the many-body 
problem is in practice never available and approximations are needed. In 
the description of CNTs we will use standard approximations which arc by 
far much more accurate and general than any approximation introduced in 
phenomenological descriptions based on model systems. 

In part I the general many-body problem is introduced. In particular 
Density-Functional Theory (DFT) and Many-Body Perturbation Theory 
(MBPT) are described according to our needs for the forthcoming parts. In 
part II the problem of double excitations is presented together with exper- 
imental evidence and the state of the art. In this part we will propose a 
new approximation which could be used in the standard approach for the 
description of absorption spectra in both the MBPT and DFT framework. 
This approximation is able to describe double excitations. Finally in part 
III the general problem of CNTs in magnetic fields will be considered. After 
a brief overview on the main experimental evidence of the Aharonov-Bohm 
effect, the Zone Folding Approach (and the Tight Binding model) will be 



viii 



Introduction 



introduced. Then we will describe how magnetic field effect are included 
in our ab initio approach and finally we will compare our results with the 
predictions of the models. 



Notation and conventions 



This is a brief overview of the conventions used to express operators and 
their related physical quantities. The same conventions are introduced in 
Ch. m Atomic units are used in Part I of the thesis, while in Part II the 
international system (SI) of units is used. 

The one body operators are written in second quantization according to 
the following expression: 

i = ^ / rf^Xid^X2dtidi2 ^^Tl<T2(xl^l,X2^2)^/'CTl(xl,^l)'^/'CT2(x2,^2)■ 



IIere (T is a spin variable, while x and t are space and time variables re- 
spectively and Acrjcr2 (xi^ii ^2^2) is the kernel of the operator. A compact 
notation is also often used 



The compact notation will be preferred wherever it will not be source 
of confusion. In this notation repeated primed variables are supposed to be 
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integrated, i.e. 

2) = iG{2', 1)W{1, 2')r*(2, 2'; 1') , 

means 

S*(l,2) = J dl'd2' iG(2',l)iy(l,2')r*(2,2';l'). 

The symbol / dV stands for J du'dt' . However, when only partial inte- 
gration will be performed (i.e. only on space, time or spin variables), the 
sums / integrals will be written explicitly. The notation G(l,2+) will be 
used for 

lim G(xi, X2, t + e), 

and the notation (A) for the expectation value on the ground state of an 
operator: 

(^oiii^o). (1) 

For the coulomb interaction we will use 

|xi - X2I 

for the one particle part of the Hamiltonian 

i?o = -iv2 + y,(x). 

A functional will be expressed with the following notation E\f^ which 
means that the energy is a functional of the density /5(x, t). 

Finally we list here some (but not all) of the symbols used in the thesis: 

• B is the magnetic field, 

• A is the vector potential, 

• $ is the magnetic flux 



xi 



• H is the magnetic induction field, 

• E is the electric field, 

• D is the electric displacement field, 

• P is the total polarization, 

• M is the total magnetization, 

• j is the current-density, with j'^P^ the paramagnetic and j^^^ the dia- 
magnetic component 

• p is the density 

• e is the dielectric constant, 

• q: is the polarizability, 

• X is the response function, xo the independent-particle one, and xks 
the Kohn-Sham one, 

• i is a four-point response function, while L is four-point in space and 
two-point in time, 

• T is the time-ordering operator while, 

• T and T[p\ arc the kinetic energy operator and the kinetic energy of 
the system respectively, 

• S is the self-energy, 

• YiH — ^ + vh, that is the sum of the self-energy and the Hartree 
potential (vh), 

• S* is the reduced self-energy, 

• Ss is a static self-energy while is a dynamical self-energy. 
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• n is the polarizability, which correspond to the density-density re- 
sponse functioijfl, 

• n* is the reduced polarizability, 

• r is the vertex function and T* is the reduced vertex function; 

• is the screened coulomb interaction, 

• g is the independent -particle Green's function, gn the Hartree Green's 
function, and G the many-body Green's function, 

• Z is the renormalization factor, 

• Exc is the exchange-correlation energy, v^c the exchange-correlation 
potential and fxc the exchange-correlation kernel, 

• S is the kernel of the Bethe-Salpeter equation (BSE) while K the 
kernel of a generalized equation (including the dynamical BSE) for the 
four-point response function, 

• /i is the magnetic susceptibility or the chemical potential, will be clear 
from the context 

• ^' the many-body wave-function and ^'o the many-body ground state, 

the Kohn Sham ground state. 

• (i>o will be used for the non-interacting many-body ground state or for 
the magnetic flux quantum $0 = 7, with h the Planck's constant and 
e the electron charge, 

• when explicitly specified 'i/s will be the wave-function in the Schrodinger's 
picture, in the Heisenberg's picture, and in the interaction pic- 
ture, 

•'The two quantities differ only because 11 is T-ordered, while Xpp is a retarded quantity 



SA stands for a variation of the quantity A, while 5{1, 2) is the Dirac's 
delta function, 
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Theoretical background 



Chapter 1 

Many— Body Systems 



Any known many-body system is constituted of interacting particles and, 
at least in principle, any physical aspect can be understood describing their 
dynamics. Elementary particles are in general ruled by the equations of 
quantum mechanics and special relativity (or general relativity) and four 
possible kind of interactions arc known to exist: the Electromagnetic inter- 
action, the Nuclear Weak interaction, the Nuclear Strong interaction and 
the Gravitational interaction. For this reason even the description of a sin- 
gle atom, where in principle all forces have to be taken into account, appear 
an almost impossible problem. Moreover any macroscopic body is consti- 
tuted of an enormous number of interacting particles, for a reference the 
Avogadro's number; Na = 6.0221415 x 10^'^ [particles/moles]. So any 
attempt to solve the many-body problems seems doomed to fail. Despite 
this discouraging scenario there are two factors which in fact make possible 
to tackle the many-body problem from a microscopic point of view. First 
the majority of the processes which happen in everyday life involve a thin 
energy window such that only the electromagnetic interaction plays a role 
and, quite often, only the dynamics of the electrons needs to be described. 
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Secondly the majority of the macroscopic objects are constituted by funda- 
mental building blocks which almost completely determine their properties: 
these are the molecules which constitutes the gases and the liquids and the 
unitary cells which are repeated an infinite number of times in many solid 
systemf0. 

A crucial role is, then, played by the equation which describes the dy- 
namics of few interacting electrons immersed in the Coulomb potential of the 
nuclei that, in a first approximation can be considered frozen in their instan- 
taneous positions: the Scrodinger Equation (SE) in the Born-Oppenheimer 
(BO) approximation. Possibly the electrons can interact with external fields 
which can be used either to explore or to tune the properties of the materi- 
als. The first part of this thesis will be dedicated to the description of light 
absorption experiments where an external light source is used to investigate 
the optical properties of a many-electrons system. The second part will be 
focused on the study of magnetic field effect on Carbon Nano- Tubes (CNTs) 
and how the external field can be used to tune the electronic properties of 
the CNTs. 

The SE can be obtained applying an Hamiltonian operator to the elec- 
tronic wave-function. The operator can be divided into two terms JJ™* -|- 
fjext 'YliQ first term describes the electrons-nuclei interaction, while the 
second term is due to the presence if an external perturbation that, in this 
work, is the electromagnetic field. The SE in atomic units reads: 



dt J 2 



, 2 

, ext 



(1-1) 



yext j^ext accouuts for thc external potentials and the term 'B^^^tr 
accounts for the interaction of the spin with a possible applied external 
field. <T is a vector constituted by the Pauli matrices, ^ is the electronic 



^For gases and liquids only the properties related to the electronic dynamics can be 
studied from a microscopic point of view. Other properties require a statistical description 
of the system. 
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wave-function; w is the Coulomb interaction while Vj accounts for the ionic 
potential. Finally the V operator takes into account the kinetic energy of 
the electrons. 

It is important to observe that some of the terms discarded in Eq. 
does not have in practice a role in the physical process we are interested 
in. The reason is that the energy scales involved are so different that the 
corresponding dynamics can be neglected. This is the case of the nuclear 
forces and of the gravitational force. Some other terms instead are discarded 
as they are usually very small although they could be needed for the descrip- 
tion of some physical phenomena. This is the case of relativistic correctionj^ 
which are needed in the description of materials composed of heavy nuclei; 
the dynamics of the nuclei, which is neglected in the BO approximation, 
is relevant for example for the description of the superconducting phase of 
some materials. Many development of the state-of-the-art are devoted to 
overcome such approximations. 

On the other hand the idea of describing solid state devices as an infinite 
repetition of the same fundamental building blocks can be applied to the de- 
scription of specific kind of materials only. As predicted by Feynman times 
ago, objects at the nanometric scale can display very peculiar properties 
which are completely different from the case of molecules or of bulk systems. 
In these direction the state-of-the-art tools need to be pushed beyond their 
present limits. From one side the common approximations involved, which 
are often based on physical intuition, cease to be valid and new approxima- 
tions are needed. From the other side the lack of a repetitive structure calls 
for the need of instruments able to describe systems with hundred, thousand 
and even more interacting electron. Only thanks to the very recent increase 

■^The terms due to the presence of a magnetic field included in Eq. I lf .f I I are already 
relativistic corrections to H"^^. Here we consider these terms as in this thesis we are 
interested in the description of the Aharonov-Bohm effect in CNTs. The relativistic 
corrections to the H^^^ part of the Hamiltonian, that is the magnetic field generated by 
the electronic current (and spin) and the magnetic field due to the nuclei (the spin orbit 
interaction mainly) are neglected here. 
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of computational power and, in the same time, to recent developments of 
the techniques, it has become possible to tackle, at least in some cases, the 
description of such system "ab-initio" and so to test the prediction of more 
simple theoretical models. 

1.1 Looking for the ground— state 

The solution to Eq. p.ip is the main goal of the many-body physics. The 
first objective is to solve such equation with H'^^^ = 0. That is to find out 
the ground-state of the system. 

An exact analytical solution can be obtained only in oversimplified sys- 
tems such as the hydrogen atom. For any realistic system this is far beyond 
our possibilities. A computational solution can be obtained at the price of a 
computational time, that grows exponentially with the number of electrons, 
only for very small systems. This is why approximations are needed. The 
problematic part of the Hamiltonian is the interaction term w(|ri — r2|) 
for which different possible strategies are available. In this work we will 
tackle the problem using Many-Body Perturbations Theory (MBPT) and 
the Density-Functional Theory (DFT). Both methods start from the con- 
sideration that the many-body wave-function contains much more informa- 
tion than really needed. So instead of the exact wave-function, which is a 
function of 3N variables, where N is the number of electrons, one can look 
for simpler quantities which contains only the informations needed to give 
a quantitative description of experiments. In quantum mechanics any mea- 
surable quantity is related to an Hermitian operator. In particular one body 
operators can be written in second quantization as: 

A=Y^ y d^xid3x2dtidt2 A^ia2(xiti,X2i2)V'ii(xi,ti)?A<T2(x2,t2), (1-2) 
We introduce here the notations 1 = (1,^) = (^^^^cr) which will be used 
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Figure 1.1: Total energy minimization: by looking at the minimum of the total 
energy is possible to obtain the ground-state properties of the system. Here, as 
an example, the total energy of Silicon Bulk as a function of the lattice parameter 
is shown. From ref. [4]. 



from now on. So for example Eq. (|1.2|) can be written in the two forms: 

A = y"dld2 A(1,2)V'^(1)V>(2) , (1.3) 
A = J dld2 A{lh,2t2)^ljH'i-h)ii2t2). (1.4) 

Within this notation we finally define 1+ = /img_j.Q+ (1, ti + e). 

Eq. (|1.2[) defines a "one body operator", that, in the jargon of the sec- 
ond quantization, acts only on one particle states. The field operator 'ijj{2) 
destroys a particle at the space-time coordinate 2 and his Hermitian coun- 
terpart ^/'^(l) creates a particle at 1. If we set X2 = xi, t2 = ti or (J2 = (Ji 
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the operator is said to be local in space, time or spin respectively. 

As observables are related to the expectation value (A) , it is then natural 
to construct a theory which is able to evaluate, instead of the many-body 
wave-function, the expectation value of any given one body operator. This 
can be done introducing the one body Green's function (GF) 

G(l, 2) = -z($o|T [V'(l)^^(2)] |<i>o>. (1.5) 

<J>o is the interacting ground-state many-body wave-function of the system 
and T is the time ordering operatoiEI- 

A particularly meaningful physical quantity is the total energy of the sys- 
tem, as its knowledge enables to obtain many informations on the system. 
An example is shown in Fig. (jl.ip . where the total energy is used to de- 
termine the equilibrium lattice parameter of silicon bulk. Unfortunately the 
total energy operators involves a two body operator, the coulomb interaction 

w{l,2) = - 1 :6{t2-h)S„,,, (1.6) 

|xi - X2I 

whose average cannot be easily obtained in terms of the GFs. Nevertheless 
the Galitiskii-Migdal equation [S] ensures that the total energy of any system 
can be expressed in terms of the one particle GF. Wc do not give the proof 
here but we just observe that, from the SE it's possible to obtain the identity 

{^^^-Ho) = {w) (1.7) 

where we have introduced Hq{1) for the one particle part of the Hamiltonian. 
All the operators on the right hand side of Eq. (|1.7p are one particle operators 



'^The T operator order any couple of operators. T[A{ti) B {t2)] = 9{ti -t2)A{ti)B{t2)± 
6{t2 — ti)B{t2)A{ti) where the sign is -I- or — according to whether the operators are 
bosonic or fermionic. 
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and so the average can be expressed in terms of tlic one particle GF: 



In the next chapter we will briefly review how the GF can be used in 
practice to get approximate expressions for several obscrvables. 

1.2 Perturbing the ground— state 

The solution of the many-body problem with i/"^-^* = provides informations 
on the system such as its mechanical stability, the electric character (metal / 
semiconductor / insulator), etc. As long as the external Hamiltonian is zero 
the body will remain in is ground-state and will not do anything special: 
a very boring condition. What is much more interesting, and at the same 
time closer to what happens in an experimental situation, is to see how a 
system react to an external perturbation. A system is never isolated and 
always interact with the enviroment. Experimentalists have to probe with 
an external field to see if a system is a conductor or an insulator, have 
to shine light on it to see which frequency the system absorbs (roughly 
speaking to see the color of the system) or they have to try to break it too 
see how much it is mechanically resistant. As a first guess one could imagine 
that knowing the ground-state we can understand how it will react to an 
external perturbation, at list to first order. This is not always the case, the 
interaction is still there and things are much more complicated (and much 
more exciting!). 

In this work we will only deal with electromagnetic external perturbations 
that we will artificially distinguish in two classes: the static perturbations 
and the time dependent ones. Static perturbations can be embodied in the 
ground-state Hamiltonian and the same techniques we will describe for the 
j^ext _ Q ^g^gg .^jjj apply. This is what we will to study on CNTs immersed 




(1.8) 




occupied states 



Figure 1.2: Perturbing tlie ground-state: a scliematic illustration of an external 
electromagnetic field which perturbs the ground-state of a semi-conductor. Elec- 
trons are excited to valence states. This is an independent particle representation. 
From ref. [5]. 



in external magnetic fields. On the other hand time dependent perturba- 
tions are more problematic to deal because, when the Hamiltonian is time 
dependent, it's not possible to speak of a ground-state of the system. New 
techniques are needed to follow the time evolution of the system; anyway 
often experimentalist just uses small external perturbations to probe the 
system. In these situations linear response calculations give us all the in- 
formation wc need to describe the experiment. For this reason other key 
quantities we will encounter are the response functions. 



The equation of motion in the presence of an electromagnetic field can be 
expressed in terms of the density, the magnetization density and the current 



1.2 Perturbing the ground state 
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density 

0"1 ■(72 
(71 ■(72 
<7\.(72 

where 

j(x,t) = ^ y'dV(^J(x-x')|^ + |^5(x-x')) +A(x,i). (1.10) 

The expression of the external Hamiltonian is then 

i?"^*(x, t) = V^'--^*(x, t)p(x, t) - A^^*(x, t)j(x, + B'=^*(x, t)m(x, t). (1.11) 

Note that Eq. (|l.lip defines a local one particle operator. Its effect on 
the electronic wave-function can be highly non linear. This is the case, for 
example, of multi-photon excitations where the one particle operator iJ^^^* 
act n-times on the ground-state wave-function, thus inducing, at least, an 
n-particles effect. In the linear regime however, only terms to first order in 
j^ext jjgg(j considered; one could be tempted to say that, then, only 

one particle process are involved. We will see that this is not the case due to 
the effect of the particle-particle interactions which induce many-particles 
process, such as double excitations. There is then an overlap between linear 
and non linear processes. This happens because the electrons react to the 
total electromagnetic field, that we have artificially divided into and w 
(and Vf); even when is small the changes induced in w could be big 
enough to induce effects beyond the linear order. We will describe double 
excitations more in details in the Part II of the present work. 

As our perturbation couples to the density, the current and the spin it is 
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convenient to look at their variation induced by weak external fieldo 

Sp = Xpp ^y^^* + Xpj <5A^^* + Xpm <5B«"* , (1.12) 
- Xjp SV^""' + Xii + Xjm , (1.13) 

5m = Xmp SV^'-'' + Xmj <^A^-* + Xmm ^B^^* . (1.14) 

5p, (5j and 5m represent the variation of the density, the magnetization- 
density and the current-density. For example 5p ~ {p{x,t)) — pQ{x,t), how- 
ever we drop the spatial and the time dependence in order to have a compact 
notation and focus on the relation between the physical quantities and ex- 
ternal perturbations. Xpp-> Xjj ^^'^ Xmm are respectively the density-density, 
current -current and spin-spin response functions. Similarly Xpjy Xpm and 
Xjm are "mixed" response functions. 5V, 5A (5B represents small external 
perturbations. Eq. (jl.l2p can be seen as the non relativistic limit of the 
general equation which couples the four potentials (V, A) to the four cur- 
rent (p, j). The magnetization enter as a result of the reduction, in the non 
relativistic limit, of the Dirac equation to the Schrodinger equatior@. 

In this thesis we will study the response to electric fields, neglecting 
much smaller dynamic magnetic field. We will work in the coulomb gauge 
to remove the coupling of the external potentials with the currents. We are 
then reduced to consider only the density response function in presence of a 
scalar potential. 

5p{l) = xil,l')5V{l'). (1.15) 

Here repeated primed variables arc integrated. From now on we will keep this 
notation. The density response function is obtained by summing over the 

*Here the expression SA = XA B^B stands for (5j4(x,t) = XA s (''i x', t')(5B(x', t'), 
where primed variables are integrated. 

^In the relativistic formulation the spin is naturally included in the four-current density. 
In particular the spin is alwrays related to a spatial current which, inserted in the Maxwell 
equations, generates the magnetic field usually associated with a magnetic spin moment. 
This magnetic field is very small and, consequently, it is generally neglected in the non— 
reativistic limit. 



1.3 The macroscopic and the microscopic world 
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spin variables: Sp{x,t) = J2„Sp{l); Xpp{^i,ti]^2,t2) = E^i<t2^(^'2) ^^^^ 
(5V^'"-^*(x,t) = 1/2E<T^(1)- The response function x(l,2) can be related 
to the Green's functions observing that it involves the expectation value of 
four field operator, two coming from the density and two from the external 
potential. In the next chapter we will derive this relation. As we can argue 
from the number of operators the relation will involve a four point, that is 
a two particle, GF. 



1.3 The macroscopic and the microscopic world 

While for the ground-state we have a receipt to directly obtain the expec- 
tation value of any physical quantity, for the linear response regime there is 
no a direct equation which relates the response function to any measurable 
quantity. 

For example in an absorption experiment, the macroscopic electromag- 
netic fields are measured. Therefore we need relations between the micro- 
scopic and macroscopic quantities. We star from the macroscopic Maxwell 
equations 



dB 

VxE = - — , VB = 0, 
ot 

WD = p™"", V X H = + 

ot 



(1.16) 



Eqs. ()1.16p together with the relations 

D = eoE + P H = ^o(B + M) (1.17) 

completely define the macroscopic fields acting in a medium. Here (E,B) 
are the electric and magnetic field, (D,H) arc the electric displacement 
and magnetizing field, (P, M) arc the polarization and magnetization of the 
system. In our work we will only consider linear systems with /i ~ /iQ. This 
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%(r",r') 

^6pind(r") 



Polarizable medium 



Figure 1.3: Polarization of a bulk system. An external charge placed in ro 
perturbes the system throught an external potential. This induce a variation of 
the density which creates an induced Electric field. From ref. [3]. 



means that: 



P = XeE, E = e-iD, 

M XmB, B = /ioH. 



(1.18) 



Finally wc consider isotropie systems where both Xe, Xm and also are 
diagonal tensors and behave as simple costants. 

D does not depend on the internal charges of the system and it can 
be identified with the external electric field, D ~ eoE"^^*. On the other 
hand E is the macroscopic total electric field, obtained as an average of the 
microscopic field over a region of space V which is big enough to smooth the 
strong oscillations of the atomic fields but small enough so that A ^ V^^^, 
where for periodic systems V correspond to the volume of the unit celjf]. 



^For absorption experiments typical energy range are around l/20eF so that A ranges 
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The macroscopic field is tlien defined as a function of tlie position of the 
unit-cell, let's call it R, and if the unit cell is small enough we can consider 
it a function of a continuous variable E(R) ~ E(x). 

In the long wave-length regime we can just take the microscopic version 
of eq. (|1.18p and average over the volume V. As the external field is slowly 
varying over the volume we obtain 

(E„,iero> = (e-Lo)D/eo (1.19) 
E = £-/ D/eo (1.20) 

(1.21) 

So all we need to simulate realistic experiments is to calculate a microscopic 
dielectric function e~Jcro which will be averaged over the volume V. To be 
precise in absorption experiments the ratio between the total field and the 
applied external field is measured, that is 

EM = — = 1 • (1.22) 

\ micro I 

The quantity e^/ is instead measured in Electron Energy Loss Spectroscopy 
(EELS) experiments where the energy loss by a fast electron, approx- 
imated as a classical particle, while traveling through a solid is measured. 
The difference between the two experiments is related to the long range term 
of the interaction as explained in App. IA.2I 

For linear isotropic materials the dielectric function is a number and, 
therfore, it can be easely manipulated in the Maxwell equations as well as in 
the equation for the electromagnetic potentials. In particular, in the coulomb 
gauge E = — Vy, and the microspopic screening can be found relating the 
external and the total potential: 1^*°* = ^-lyext ^ .^jjj define Smicro 
in the next chapter and we will see that it is related to response function 



from lfj,m to 50nm 
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through the relation 

C.L-o = l + ^X- (1.23) 

In the case of isolated systems it is not possible to define a meaningful 
and finite volume V. Strictly speaking — > oo and all averaged quantities 
goes to zero. Consequently, instead of the dielectric function, we introduce 
the polarizability a which relates the change of the dipole moment to the 
external field EU. 

6p^a— (1.24) 



Here the long wave-length limit appear as a dipole expansion for the 
external applied field. As fc = ^ << 1 the external field and the external 
potential can be written as: 



E'=^*(x, t) = Eoe"'''""* ~ Eoe"'"*, (1.25) 



from which 



l^'=^*(x,t) ~ -Eoxe"^"* (1.26) 



Observing now that Sp = —{xSp) we can construct the microscopic equiva- 
lent of cq. (|1.24p starting from equation Sp = xSV'^'^*. 

(x<5p) - (xxx)Eo (1.27) 

with a = (x^x) the polarizability tensor. Although isolated systems are al- 
most never isotropic, experimentally the absorption spectrum of molecules is 
often obtained in the gas phas(j§. The gas by itself is isotropic and the exper- 



^ Microscopic fields do not appear in the macroscopic equations as their macroscopic 
average goes to zero. However they are still present. Here for example microscopic fields 
modify the dipole. 

*Many measurements for molecules are in condensed phase, either in solution or in 
molecular crystals. For these statistical analysis is carried out to relate molecular response 
properties to macroscopic response properties. Solution and solid effects are more and 
more frequently included in calculations of molecular absorption spectra. 
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iment can be interpreted by using the spatial average of the polarizability: 
Tr{a) = {axx + ayy + a^z)/3. 
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Chapter 2 

Green's Function approach 



In the previous chapter we have introduced the one particle Green's Function 
(GF) G{1,2) as the natural quantity to describe in order to evaluate the 
ground state properties of the system. In contrast to the Many-Body wave- 
function, which is an 3-/V + 1 variables function, with A'^ the number of 
electrons, G is a function of only two space-time coordinates, i.e. 8 variables. 
Nevertheless the GF of a fully interacting many-body system is as complicate 
to calculate as the ground state wave-function. This can be easily seen 
writing down the equation of motion for G, starting from the relation 



dA{t) 



(2.1) 



dt 



where A{t) is any operator. In the case of the field operators '0(1), V'^(l) 
the Eq. (|2.1|) yields an equation of motion for the GF 



i- - Ho{l) G(l, 2) = <5(1, 2) - tU{l, 3)G(1, 3; 2, 3+). 



(2.2) 
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Eq. (|2.2p shows that the one particle GF depends on tlie two particles GF. 
In the same way the two particles GF introduces the three particles and so 
on. This defines an hierarchy of equations that cannot be closed exactly. 
Approximations are required. 

There is an advantage in computing the GF instead of the many-body 
wave-function: the GF is the minimal object needed to compute the expecta- 
tion value of any one particle operator, while the many-bodywave-function 
contains much more informations which often are not ncedccl^. In this chap- 
ter we will introduce the techniques to evaluate G(l, 2) and we will illustrate 
how the response function, the key quantity of linear response theory, can 
be written in terms of the GF. 

2.1 The time— evolution operator and the role 
of the interaction 

The GF approach allows, through the second quantization formalism, to 
isolate the complicated part of the Hamiltonian, the interaction. This can 
be done by writing the field operators in the "interaction picture" . In the 
Schrodingcr picture the SE reads: 

^^^\M,sit))=msit)). (2.3) 

The wave-functions are the key quantities of the theory and it is possible to 
define a time evolution operator U{t) ~ e~*^*. In second quantization the 
key quantities are the operators, so the idea of the Heisenberg picture is to 
describe the system evolution in terms of the operators. The wave-functions 



^Only for small molecules wave-function methods remain dominant as they give more 
accurate answers than Green's function methods. However for infinite system explicitly- 
correlated wave-function methods are impractical. 
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is imposed to be static: 

\^h) = e'^*|*sW), (2.4) 

and accordingly the time evolution operator correspond to the identity. By 
demanding that 

(*a|6H|*ff) = (*s|0|*s), (2.5) 
the time dependent operators are, then, constructed: 

Oh = e'"*de~'"*. (2.6) 

The interaction picture is obtained in a similar way by moving only the 
complicated part of the time evolution on the operators. So the Hamiltonian 
is split in two parts: 

H = Ho + Hi, (2.7) 
and the wave-function is defined as 

\^j{t))^e'"''*\^s{t)). (2.8) 

As a consequence the wave-functions arc solutions of the following Schrodinger 
equation 

«^|*/(i)) = e'""'Hie~'""'\<fs{t)). (2.9) 
Similarly operators are written as 

Oi = e'^°*de-'"°K (2.10) 
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In the interaction picture the time evolution of the wave-function can be 
cast in a generalized time evolution operator 



+00 „t „t 



n=0 

T 



(2.11) 



defined in such a way that 



\^i{t)) ^Ui{t,to)\^{to)). 



(2.12) 



2.2 Equilibrium properties 
The one particle propagator 

The interaction picture is a practical starting point thanks to the Gell-Mann 
and Low theorem[5]. The Hamiltonian 

H,{t) ^ Ho + e-'^'^Hi (2.13) 

is considered By using Eq. (|2.12p the time evolution of any cigenstatc can 
be written as 

mt))-M^^^^M^, (2.14) 

($f|C/e(t,io)|$}) 

where |$}) is the eigenstate at time t = t^. The Gell-Mann and Low Theorem 
states that in the limit to ^ —00 and e — > o[§ $} reduce to an eigenstate 
oi Hq + Hi. A relation among the eigenstates of Ho and the eigenstate 
of Ho + Hi is then established. Finally assuming that the not interacting 
ground state slowly evolves to the ground state of the interacting system we 



^It is important to take the limit in this order in order to obtain a meaninful result. 
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obtain [5] «The most useful result of quantum field theory» 

n J-oo 



i?i(ti0...i?i(t„'W(l)^/(2) |$o) 



($o|^|$o) 



(2.15) 



(c&olT pl(l)V^/(2)J |<f>o> 
' (*o|^|*o) ^ 



with 



4-00 



dty . . . dtn'T 



Hi{ty) . . . Hiit„,) 



(2.16) 



Eq. (|2.15p is a perturbation expansion in the interaction which involve the 
expectation values of many field operators. Any term of order n involves 
n-times the interaction operator Hi = w(l, 2). Thanks to the Wick theorem 
[S] we can express each term in the scries in terms of the not interacting GF 



«ff(l,2) 



(*o|T 



V'l(l)Vi7(2) 



l*o) 



($o|*o) 



(2.17) 



Any term in the expansion starts contains a not interacting GF which start 
from point 2, g{ ,2), and a not interacting GF which ends up at point 
1, 5(1, ), while the presence of the interaction affect the propagation in 
between. So the interacting GF can be expressed as 



G(l, 2) = g{l, 2) + g{l, 1')Sh(1', 2')5(2', 2), 



(2.18) 



where we have introduced the I]i/(1,2) self-energy function. This functions 
includes terms which are repeated an infinite number of times. By introduc- 
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ing a reduced self-energy SJj (1, 2) which satisfy the equation 

^h{1, 2) = 2) + EJ,(1, l').g(2', 2')Sif(2', 2), (2.19) 

Eq. (|2.18p then becomes 

G(l, 2) = .9(1, 2) + g{l, l')^*„{l\ 2')G{2', 2). (2.20) 

Eq. (|2.15p and Eq. (|2.20[) represent the starting point for the development 
of approximations to the many-body problem. In particular Eq. p.l5|) can 
be approximated by truncating the scries expansion to finite order. This is 
a reasonable approximation for systems where correlations effects are less 
important and the interaction can be treated as a small perturbation. For 
more correlated system, and in general for solids, the expansion in the bare 
interaction in not meaningful as screening effects are dominant. As discussed 
in the next sections in this regime the electrons are screened and the interac- 
tion is considerably weaker than in isolated systems. Eq. ()2.20p can be used, 
choosing approximations for the self-energy EJ^ , one of the most common 
being the GW approximation, which is the first order term in a possible 
expansion for the GF in powers of the screened interaction W. 

Finally inserting Eq. (|2:20|) in Eq. ([22]) the EOM for the GF reads 

G(l, 2) = S{1,2) + SJ,(1, l')G(l', 2). (2.21) 

The Hedin's equations (I): bare interaction 

To device useful approximations for the self-energy we follow the method 
developed by Lars Hedin [5]. An expression for the self-energy can be ob- 
tained by studying how the full GF react to a fictitious external potential 
S(p{l). This potential is then set to zero at the end of the derivation to 
recover the ground state GF. 



^^-Ho(l) 
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In the following we include the interaction in the first part of the Hamil- 
tonian H = Hq + w{1,2), while the potential ^i(l) = ip'' {^)'4'{^)^v{^) is 
included in the S operator: 



iG^(l,2) 



*o|T 



5^t(i)^(2) 



(*o|5|*o> 

The variation of the GF with respect to the potential Sip reads 



-G^(l, 3; 2, 3+) + G^(l, 2)G^(3, 3+). 



(2.22) 



(2.23) 



Eq. ()2.23p can be used in the Eq. p.2p . to express the two particle GF in 
terms of G(l,2): 



^^-i^o(l) 



G(l, 2) + i w{l, l')G(l', l'+)G(l, 2) 



Using the equality 

5G(1,2) 



we define the vertex function 



r(l,2;3) 



5ip{3) ■ 



(2.25) 



(2.26) 



By using Eq. (|2.2ip and Eq. ()2.26p we obtain the following expression for 
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the self-energy: 

2) = -tS{l, 2) w{l, 1')G{1', !'+) + 

— (1,2')G(1,1')^^^ (2.27) 
= S{1, 2)vHil) + t w{l, 3')G(1, l')r(l', 2; 3'). 

Using tiie Dyson equation for 0(1,2) to compute the functional derivative 
of we get 

G(l, 2) = gnil, 2) + .9^(1, l')S*(l', 2')G(2', 2), (2.28) 

I]*(l, 2) = *G(2', l)w{l, 2')r(2, 2'; 1'), (2.29) 
r(l,2;3) -<5(1,2)5(1,3)+ 

The Hartree potential is included in the definition of the bare GF 

2) = ,9-1(1, 2) + z;^ (1)<5(1, 2). (2.31) 



The Hedin equations (II): screened interaction 

Lars Hedin introduced a fundamental breakthrough by realizing that, instead 
of the bare Coulomb interaction ti;(l,2), one should consider a perturbative 
expansion in the screened potential W{1, 2). In many-body systems, indeed, 
the interaction is always screened. Formally we introduce the potential 

SV{1) ^ 5^{1) + 5vh{1) (2.32) 

defined as the sum of the fictitious potential ip and its classical screening. 
This is the classical screened potential because it is obtained by neglecting 
the changes in the quantistic self-energy. Considering linear and isotropic 
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systems the microscopic dielectric function can be defined as the relation 
between the total potential and the external potential (see Ch. ([T]), Sec. 
(113): 

^-'(1,2)= ^= <5(1,2)+ u;(f,f')n(f',2), (2.33) 
where 11(1, 2) is the response function 

n<— (-) 

The screened interaction is then 

iy(l,2) = 6-i(l,l')^l',2) 

(2.35) 

= w(l,2)+ u;(l,l')n(l',2')u.(2',2). 

The polarization function can be interpreted as the self-energy of the screened 
interaction W and, similarly to the case of E, it can be reduced 

^ ' ^ '5^(2') <5^(2) (2.36) 

= n*(i, 2) + n*(i, 1 V(i', 2')n(2', 2) , 

to obtain 

M^(l,2) = u)(l,2) + ?i;(l,l')n*(l',2')W^(2',2) (2.37) 
The equation for the reduced polarization can be obtained From Eq. (|2.25p 



n*(l, 2) = -i G{1', 1)G(1, 2')r*(l', 2'; 2). (2.38) 
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G(l,2) 5(1,2) .9(1,10 S^2'')G(2', 2) 

< = < + — < — Csf^)=<= 




Figure 2.1: Diagrammatic representation of tlie Hedin Equations. 



where we have defined a reduced vertex function 
G-i(l,2) 



r*(l,2;3) 



5V{3) 
(5(1,2)5(1,3) ' 



SVi^) (2.39) 
= <5(1,2)<5(1,3) + 

g^G(l',3')G(4',2')r'^(3',4';3). 

The self-energy can be expressed in terms of the new quantities in order to 
obtain a closed set of equations, the Hedin equations: 
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G(l, 2) = .9h(1, 2) + .gff (1, 2')G(2', 2), (2.40) 

W{1, 2) = 2) + mil, l')n*(l', 2')W^(2', 2), (2.41) 

n*(l, 2) = -i G(l, l')G(2', l)r*(l', 2'; 2), (2.42) 

S*(l, 2) = iG(2', l)W{l, 2')r*(2, 2'; l'), (2.43) 
r*(l,2;3) = <5(l,2)J(l,3)+ 

^^±lJ-^G{2', 3')G(4', l')r*(3', 4'; 3). (2.44) 

This set of equations is exact. Indeed their solution is not easier than the so- 
lution of the many-body problem. Nevertheless the Hedin's equation offer a 
convenient starting point to develop approximations. One of the most com- 
mon approximation to this set of equation is to set r*(l; 2, 3) = (5(1, 2)5(1, 3), 
that is to drop the complicate functional derivative ^(^^- This approxi- 
mation is known as GW approximation. The GW approximation is simply 
the first order approximation in the expansion of the self-energy in powers 
of the screened interaction [5] and goes one step beyond the Hartree and the 
Hartree Fock approximations. The infinite order resummation of Feynman 
diagrams is already included in S which is then iterated to construct the 
GF. Indeed to Dyson equation are solved, one for w(l,2) and another for 
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G(l,2). 



2.3 The quasiparticle concept 

Starting from Eq. p.2ip it is possible to introduce the quasiparticle (QP) 
concept, which was first introduced by Landau in the description of Fermi 
liquids^. First we Fourier transform eq. (|2.2ip 

{lu - Ho{l) - vh{1)) G(1, 2; w) - l';c.)G(l', 2;lu) = S{1, 2), (2.45) 

then perform an analytic continuation of the frequency variable to the com- 
plex plane, a; — z. A formal solution to the equation can be obtained by 
using a biorthonormal representation of the GF 

The right and left wave-functions '^\{l,z) and ^\{l,z) satisfy the equations 

(Ho(l) + vh{1)) *a(1, z) + ^"{1, 1'; z)-^y,{l, z) = SaW^a^, 
(i/o(l) + t'H(l)) *a(1, ^) + 2-^^(1, 1'; 2)I'a(1, z) = Sa(z)I'a(1, z). 

The QP concept can be introduced assuming that the dominant contri- 
bution to the GF comes from the fixed points of E\{z): 

Ef'' = Ex{Ef''). (2.48) 



The QP equation is then defined as 

(iJo(l) + vh{1)) *.(1) + S*(l, 1'; = ^^^''^.(l)- (2.49) 
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The renormalization factor 

In practice the solution of eq. p.49p is usually obtained as a first order 
correction to some mean field theory, like the Hartree approximation, that 
is eq. (|2.49p with S = 0. The solution of Eq. (|2.48p is the approximated by 



(2.50) 



Eq. (|2.50p is solved by linearizing the frequency dependence of the Self- 
Energy around the Hartree poles en 



(E(£;0^)) c (S(6^)) 



(2.51) 



Defining the renormalization factor 



,9S(a;) 



(2.52) 



we obtain 



QP ~ 



ef +Z{j:{Ef'')) 



(2.53) 



In general the number of GF's poles E^ is larger then the number of 
IP states. Indeed only if the self-energy frequency dependence is linearized 
according to Eq. (|2.5ip the GF poles coincide with the QP states. Such an 
approximation fails in the description of satellites, as the one shown in Fig. 
(|2.3p . which appear as extra poles in the spectral function A{u}) defined as 

{A{uj)) = -signiji - uj){Img[G{uj)]) 

TT 

^ 1 \Img[i:M\\ ^^-^^^ 

ti{lo- e, - Re[T.,{Lu)]f + (/m[E,(t^)])2 ' 



In the second line the non diagonal part of the self-energy is neglected. 
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Figure 2.3: Quasi Particle representation of the Green's Function. The pole 
of the not-interacting GF is shifted and broadened. QP are connected to the 
coherent part of E. Additional peaks in the GF can appear due to the incoherent 
contributions (plasmons, polarons, resonances,...). These additional peaks cannot 
be interpreted in a IP theory. From ref. [6] 

The extra poles, in the quantum chemistry language, can be described as 
"multiple excitations" [10] . Eq. p.50p is a good approximation only if the 
starting wave-functions are closed enough to the QP wave-functions. This 
is usually not the case for the Hartree theory and for this reason in practice 
QP corrections are calculated on top of Density-Functional Theory (DPT) 
calculations. We will introduce DPT in the next chapter. 

We will find a similar distinction between QPs and satellites (extra-poles) 
in the case of neutral-excitations. In Sec. (|2.4p the Bethe Salpeter Equa- 
tion (BSE), a Dyson-like equation to compute neutral excitations, will be 
introduced. Like the self-energy, S, the kernel of the BSE, S, will be able 
to produce extra poles; however while E is usually almost diagonal, S is not. 
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-tU|5 



Figure 2.4: Lehmann representation, on tiie left, compared with quasiparticle 
representation, on the right. From ref. 

This means that IP poles of the GF are shifted to QP and all incoherent 
peaks are easily identified as satellites, while in the BSE IP electrons-holes 
transition are strongly mixed and the resulting spectrum is usually broaden, 
making difficult to identify possible satellites. This argument will be subject 
of the Part. II of the present thesis on double excitations. 



The Lehmann representation 

The biorthonormal representation of the GF in the complex plane is linked 
to the Lehmann representation of the GF on the real frequencies axis. This 
can be obtained using the full many-body wave-function and the identity 




(2.55) 
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with the exact many-body eigenstate. By inserting Eq. ()2.55p in tlie 
GF definition, the Fourier transform of G yields 



U) — LOj + iri 

U! — U!j — irj 



(2.56) 



Eq. (j2.56p defined the Lehmann representation. This is an exact representa- 
tion expressed on the real axis instead of the complex plane. The Lehmann 
and the biorthonormal representations are closely linked, as exemplified in 
Fig. (|2.4p . The concept represented in the figure emerge if we consider the 
thermodynamical limit of the Lehmann representation 

G(l,2;.) ^ v/.(l)/}(2)^_g.(l)g}(2) 

uj — ujj + ir] uj + (jjj — irj 

(2.58) 

UJ ~ X 

where /./(I) = (1'oK"t(l)|^,;), 5,7(1) = o|^(2)|*j) and 
M^) = E (^fjWrAmf^-Ej) + gj{l)g:j{2)e{Ej-fi)y{Ej-x) (2.59) 

is another expression for the spectral function previously defined. Eq. (j2.57p 
is an integral representation of the GF which can be analytically continued: 
w — > 2. In complex analysis a branch cut, here the series of poles in the 
Lehmann representation of the GF in the thermodynamical limit, can be 
expressed as the integral of a complex pole, here the poles of the spectral 
function. 



In this respect the QP concept well describes the power of MBPT. In- 
deed the QP poles are not eigenstate of the Hamiltonian being formed by a 
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macroscopically large number of almost degenerate eigenstates of H. QPs 
have a finite lifetime. 

An infinite series of closely lying poles which appear as a branch-cut in 
the Lchmann representation is replaced with a single complex pole in the QP 
picture. Nevertheless the QP representation fails to describe the continuous 
part of the GF. For other details on the QP concept see App. |B] 

2.4 Neutral excitations: the Bethe— Salpeter 
equation 

Once we have computed the one particle GF we can evaluate any physical 
proprieties of the system at rest. Still we have no information about how 
the system would react if disturbed by a time dependent external potential. 
In principle it is possible to calculate the perturbed time dependent one 
particle GF, but in the case of weak perturbations we can use the linear 
response. In this approach the key quantity is the linear response function 
which can be calculated in terms of the GF. This is the subject of the present 
section. The Lehmann representation shows that the linear response function 
describes the neutral excitations of the system where two or more particles 
are involved. Therefore the two particle GF will naturally emerge as the 
starting point to compute linear response properties. 

We start from the wave-function evolution induced by a time dependent 
perturbation 




(2.60) 
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Then the expectation value of any operator can be written as 




sv 



ext 



(2.62) 



and by using the density operator Eq. (|2.6ip becomes: 



p(i) = po(i) + e{t t')sv^^\r){^m [p{i),p{i')] l*(o)>. (2.63) 



Formally this is different from the polarization function 11 that is a T-ordered 
quantity while x is a retarded quantity. However it can be proven through 
the Lehmann representation that the two are equivalent a part from a small 
shift of the poles along the imaginary axis. The T-ordered response function 
is: 



where the operator Sp = p — {p)o is used to ensure that the terms which 
involves the expectation value on the ground state are canceled, as for the 
retarded function. 

Wc have already evaluated 11 by solving the Hedin equations. However 
in contrast to the GF, the GW approximation, that is F = 1, is not a good 
approximation for the response function. The reason is that we are looking 
to different physical processes. In order to find a better approximation for 



The retarded response function is thus defined as 



(2.64) 



zn(l,2) = (*(0)|r[<5/5(l)<5p(2)]|vI/(0)) 



(2.65) 



2.4 Neutral excitations: the Bethe Salpeter equation 
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n we use the Hedin equation for the vertex 

5 

dG{l',2') 



r(l, 2; 3) = ,5(1, 2)5{1, 3) + ^^7^ [^(1, 2^(2) + 2)] 



G(2', 3')G'(4', l')r(3', 4'; 3), (2.66) 
to obtain an equation for the Polarization 

n(l, 2) = -iG(l, 2)G(2, 1) + -iG{l, l')G(2', 1) 

^^(3,^^,) [-5(1', 2')vHi2') + 2')] G(4', 5')G(6', 3')r(5', 6'; 2). (2.67) 

Eq. (|2.67p has the structure of a Dyson equation of the form, 

n = Ho + HoSn, (2.68) 

with Ho = — iGG. However the combination — zGGF factor appearing in 
r.h.s of Eq. (|2.67p is not a two point quantity. The reason is that the kernel, 

2(1, 2; 3, 4) = ^^^^ [-5(1, 2)1^^(2) + £*(!, 2)] , (2.69) 

is a four point quantity. Eq. (j2.68p can be closed using the extended space 
of two-particles GFs, defined as 

iL{l, 2; 3, 4) = -G2(l, 2; 3, 4) + G(l, 3)G(2, 4). (2.70) 

The polarization is written in terms of L as 

L(l,2;l,2) = n(l,2), (2.71) 
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Looking at the derivation of tlie Hedin equations we see that 

and using Eq. (j2.25p and the definition of the vertex function we finally 
obtain 

2; 3, 4) = -iG{l, 4)G(2, 3) - iG{l, l')G(2', 2) 

. 6{diV,2')vH{l') + ^*{l',2')) y5G(3^3) 
* 5G(3',4') ^ ^^¥'(4,40 

= Lo{l, 2; 3, 4) + Lo{l, 2; l', 2')S(1', 2'; 3', 4')L(3', 4'; 3, 4). 

(2.73) 

Eq. ()2.73|) is the Bethe Salpeter Equation (BSE). It's then possible to con- 
struct approximations to the response function, via the kernel S. The most 
common approximation is derived from the GW approximation to the self- 
energy, S = iGW, by neglecting terms second order in W. The kernel 
reads: 

S(12, 34) = J(l, 2)S{3, 4)w;(l, 3) - ,5(1, 3)5(2, A)W{1, 2) (2.74) 

The screened interaction appearing in S is usually taken as static. In the 
case of double excitations we will see that this static approximation cannot 
describe multiple neutral excitations. 



Chapter 3 



Density— Functional 
Theory 

In the previous chapter we have introduced the MBPT formahsm. The main 
advantage of this approach, compared to the solution of the full many-body 
Schrodinger equation, is that it makes possible to develop efficient approxi- 
mations. The reason is that the electron-electron interaction, which consti- 
tutes the complicated part of the many-body Hamiltonian, can be treated 
in a perturbative manner. The solution of the MBPT equations is anyway 
a demanding task and simpler approaches are desiderable. In this chapter 
we will present the Density-Functional Theory (DFT). The key quantity of 
this approach is the density of the system which is a considerably simpler 
quantity than the GF. The density completely describes a system of classi- 
cal interacting particles and is the only quantity needed within the Hartree 
approximation. To get a more realistic description, terms approximating 
the quantum mechanical effects need to be expressed as functional of the 
density. 

The first attempt in this direction goes back to the Thomas-Fermi (TF) 
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model[T21|T3] in the 1927 and its extension due to Dirac in 1928 [Tl|. However 
Tomas-Fermi-Dirac theory remained rather inaccurate for most apphcations. 
Only some time later, in 1964, the Hohenberg-Kohn (HK) theorem [T5] put 
the theory on a firm theoretical footing. The HK states that the ground- 
state of a is an exact functional of the density. The proof of the existence 
of an exact functional gave a boost in the developments of density-based 
approach. However this was mainly thanks to the surprising success of the 
Kohn-Sham (KS) scheme (1965) [16], within the Local-Density Approxima- 
tion (LDA), that DFT became popular. In this chapter we will start from 
the demonstration of the HK theorem and then we will introduce the KS 
scheme and the LDA approximation. In particular we will show why the 
LDA approximation, which had been initially designed for uniform systems, 
has been found to be accurate even for strongly inhomogeneous systems such 
as isolated atoms and molecules. 

While MBPT offers a direct recipe to evaluate all ground-state properties 
starting from the Green function, DFT gives only the total density and the 
total energy. As stated in the introduction of the present thesis, by using 
the total energy many properties of the system can be obtained but others, 
as for example the electronic gap, are not directly accessible. Anyway it 
is common practice to interpret the KS band gap and to look at the KS 
wave-functions to obtain more informations on the system. In practice DFT 
provides a "zero-order Hamiltonian" for a first understanding of the physical 
properties of the system and also for MBPT calculations. 



3.1 The Hohenberg— Kohn theorem 

The HK theorem states that: [HKl] for any system constituted by inter- 
acting particles, for a fixed the interaction, there exists a bijectivc relation 
among the external potential (up to an arbitrary additive constant), the 
ground-state many-body wave-function and the ground-state density of 
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the system, 

y-* ^ *o ^ Po ; 

[HK2] the ground-state energy and density can be determined by minimiz- 
ing a functional of the charge-density. As a consequence of the HK theorem 
the knowledge of the ground-state wave-function gives access to all physi- 
cal observables. Moreover the external potential fixes the Hamiltonian and 
knowing the Hamiltonian we can in principle access the excited states of 
the system [T7]. It is surprising that given the ground-state density we can 
access all this physical information. However DFT is in practice used only 
for ground-state properties, as there is no practical scheme able to describe 
excited stated starting from the ground-state density. Excited state proper- 
ties can be obtained from the extension to the time domain of DFT: this is 
Time-Dependent DFT (TDDFT). 

The proof of the [HKl] theorem is straightforward. First we show that 
given a ground-state wave-function there exists a unique external potential 
which determines it. The relation in the opposite direction is trivial, at 
least for systems which does not present a degenerate ground-state, as the 
Schrodinger equation has a unique solution. The second step is to show 
that given a ground-state density there exists a unique wave-function which 
determines such a density. Again the reverse relation is trivial as the density 
is an observable of the system and any observable can be obtained from the 
wave-function: 

Po = («'o|p|*o>. (3.1) 

Step one. Let's suppose that there exist two different external potentials 
which have the same ground-state wave-function: 

{f + w + V--')\^>o)^E„\'fo) 
{f + w + V''^*)\^o)=E'o\^o)- (3.2) 
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If we subtract the two equations assuming that the wave-functions are dif- 
ferent from zero only in a subspace of null measure we obtain 

F'=^*(x)-T/"=^*(x) = (£;o--B^,). (3.3) 

This means that the two potentials can differ only up to a constant which 
anyway is irrelevant as the total energy of a system is always determinated 
up to an additive constant. 

Step two. Suppose that ^'o and ^'q are the ground-state wave-function 
of two different Hamiltonians H and H' . From the definition of ground-state 
it follows that 

Eo = (^-oli^l^o) < (*ol^l*o). (3.4) 
The r.h.s. of Eq. p.4p can be written as 

{%\H' - 1/"=^* + V'^'^'l^) =Eo + J d^r p'{r) [y'=^*(r) - y"=^*(r)] , (3.5) 

which means 

Eo<E'o + J d^r p'(r) [y^"*(r) - r^"*(r)] . (3.6) 

By following the same procedure but inverting the roles of and we 
obtain 

Eo<E'g + J d^r p'{r) [y'=^*(r) - T/"=^^*(r)] . (3.7) 
If we assume po = p'q, summing up Eq. p.Op and p.7p we obtain 

Eo + E'o < Eo + E'o. (3.8) 

This is clearly impossible and so the assumption po = p'q is wrong. 

Once [HKl] has been provcr0 the functional of the [IIK2] part of the 



^Thc proof of HKl assumes v-representability of the density: i.e. that given a "reason- 
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theorem can be written as: 

E[p] = {m\H\m)- (3.9) 

This functional has a minimum when p = is the ground state density of 
the system. Moreover it can be written as E[p\ = F[p\ + J (i'^X(o(x)F'^^*(x), 
where F[p] — {^[p]\T + w|\E'[p]) is an universal functional for any many- 
electrons system. 

The basic assumption in the proof of the theorem is that the coupling 
of the system with the environment is given by a term of the form V'^'^* p. 
This is not always the case. If an external magnetic field is present then the 
terms Aj or Bcr or both must be included. In this case the modified HK 
theorem leads to Spin Density-Functional Theory (SDFT) for density and 
magnetization [T7], Current-Density-Functional Theory (CDFT) for density 
and current |181ll9j and to SCDFT when density, magnetization and current 
are considered. 



3.2 The Kohn— Sham scheme 

The HK theorem states that the density determines the ground-state of 
an interacting system. Still, it does not provide any recipe to use it. The 
prescription is introduced by the KS scheme. The idea is to use the HK 
theorem for an auxiliary system of non-interacting particles whose ground- 
state density is assumed to be the same as that of the interacting system: 
Ps = Po- Then we look for the external potential Vs[p\ and for the non- 
interacting wave-function ^'^[p] related to such density by the HK theorem. 



ably well behaved" non negative function p(r) one can always find a local external potential 
y^^*(r), so that p(r) is the ground state density of the Hamiltonian H = T + w + yext ^ 
Unfortunately this is not always the case, for a detailed discussion see Ref. 1171 . 
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We start from the energy functional of the real system 

E[p\ = (M/bllf + y-*+u;|vI/[p]) =r[p] + f/-*[p] + C/„,[p], (3.10) 

then we rewrite it, by using the relation p ■4=> starting from the non 
interacting many-body wave-function 



^ 1"^ I (3.11) 

^Up\ + u''''[p\ + Eh[p\ + e.,M- 

The non-interacting wave-function enters the kinetic energy term, while the 
external potential energy is the same as that of the interacting system. We 
have also introduced a term which describes the energy of an interacting 
system of classical particle Eh\p\ and a last term -E^cM so that E[p\ ~ Es[p\- 
This means that 

E^M = (T[p] - TM) + {Eu\p] - Eh[p\) . (3.12) 

Exc[p\ is the xc-energy. Now we look for the ground-state by minimizing the 
functional Es[p]. In a non-interacting system = ni|V'i0 and p = lipi]"^ 
so that we can minimize the energy functional with respect to ipt. This is 
a constrained minimization as we want the tpi to be orthonormal. That is, 
using the theory of Lagrangian multiplier, we have to minimize the functional 



E[p{^l;i,...,i;n)]+^(^5h^k- ^h.k J V'h(r)V'fc(r)rf^r 
Using the relation 



(3.13) 



S S S 



^Herc the notation Ilil *i | means one have to perform an "antisymmctrized" product 
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we find 



(3.15) 



that leads to 



{i + vh[p] + v"''* + i': 



'xc 



[p]) ipi = Xiipi. 



(3.16) 



Eq. p.l6p is the KS equation, v^dp] = 5Exc[p\/ 5p is the unknown xc- 
potentiajfl which encloses all the difhculties of the interacting many body 
system beyond the classical Hartree potential. Though Eq. p.l6p has to be 
solved sclf-consistently with respect to the density, the scheme offers a very 
appealing starting point as the xc-potential appears in the equations in the 
form of a local multiplicative operator. If compared to the exchange term of 
the HF approach or to the Self-Energy of the MBPT, which are non local 
operators, the advantage is evident. Unfortunately no hint to approximate 
the xc-potential or the xc-energy functionals is given. 

3.3 The local— density approximation 

One of the reasons beyond the success of DFT, despite its simplicity, is that 
the simplest approximation proposed for Exc^ the Local-Density Approx- 
imation (LDA) successfully describes both extended and isolated systems. 
Such an approximation was designed to work for systems where the den- 
sity is almost spatially uniform or slowly varying. The idea is to compute 
the exchange correlation energy of the uniform electron g function of 

the constant density, poj £'5°™c(po) and then to express the general energy 
functional of the density as 



^Here we assumed that the functional derivative exists; however this must be verified 
for any given energy functional. 




(3.17) 
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where 

4rip) = ^"yf^ (3.18) 

is the energy per electron. 

The total energy of the homogeneous electron gas can be obtained accu- 
rately from Quantum Monte-Carlo calculations. 



The success of the local— density approximation 

To explain the success of the LDA we need to introduce the exchange- 
correlation hole-density which is linked to the success of the LDA approxi- 
mation. 

First we introduce the expectation value of the pair correlation function 

,H(.,0^<«-fc:i2. (3.19) 

p(r)p(r') p(r) 

which can be used to write the energy correlation due to the interaction (w) 

as 

EM = lld'rl d^r' ^^:l^g[p](r, r'). (3.20) 

Than we make use of the adiabatic connection and of the HK theorems. We 
define a group of Hamiltonians Hx such that 

Hx^f + Y, [ d^'v vx{v)i^l{v)4,„{v) 

+ ^E/ J d'rd^r' w{r,r')^Pl{r)4>l{r')Mr)i'<r'{r'), (3.21) 

with vx such that px = Pi for each A. For A = we have the KS system. 



3.3 The local density approximation 
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-o.E -ai CLi -1.0 -oa □ cls 



(a) Exchange-correlation hole 




(b) Spherical average of the Exchange— correlation hole 



Figure 3.1: Exchange-correlation hole and its spherical average within LDA are 
compared with the exact solution at two different densities r/ao — 0.09 and r/ao = 
0.4. From ref. [20]. 
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Now we compute, using the Hellmann-Feynman theorem 



dEjX) 
dX 




(3.22) 



where gA[p](r, r') is the pair correlation function for the Hamihonian H\ and 
its ground-state. Finally we express the total energy for A = 1 as 



E{1) = E{Q) + [ dX 
Jq 



dE{X) 
dX 



d^rdh' w{r, r')p(r)p(r') / dX {gx [p] (r, r') - 1) , 



(3.23) 



where, in the last line, we have subtracted the Hartree energy from the 
expression for {'w)x- This last term is clearly the exchange correlation energy. 
Defining the exchange-correlation hole-density 



gxc{r,r') = p(r') 



dX{gx[p]ir,r') 



11 



(3.24) 



we can therefore study the behavior of the g^c to check the LDA performance. 
In particular Fig. p.ip shows how g^c behaves for the neon atom. The exact 
solution is compared with the LDA one. In the two upper frames, (a), g^c is 
plotted for a fixed r coordinate, as a function of r' along a direction parallel 
to r. In the two bottom panels (6) its spherical average is plotted. We see 
that compared to the exact xc-hole the LDA performs quite badly, while the 
gxc spherical average is well described. For an interaction which depends 
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only on the modulus of distance among the particle Exc[p\ depends only on 
the spherical average and this is one of the reasons why LDA well performs 
even for not homogeneous systems. Another reason is that the LDA xc-hole 
satisfy the sum rule 

y"dV^.,,(r,r') = -1, (3.25) 
which helps to guarantee error cancellations. 



Spin dependent local density approximation 

For spin polarized systems, that is systems whose ground-state present a 
magnetization jti^ 7^ 0, the DFT approach could be ideally used. However 
in these approach the magnetization would be an unknown functional of 
the density and the LDA approximation is not guaranteed to work. For 
these reason it is preferable to work within the SDFT formalism where the 
magnetization can be obtained directly from the KS wave-functions and its 
possible to design the Spin LDA (SLDA) as a direct extension of the LDA. 
To this end the quantity 

Pt + Pi 

which measures the polarization of the system, is defined. Then the approx- 
imation 

e.cip, = e^c(p, C = 0) + (e,e(p, C = 1) - e.c(p, C = 0)) .g(0 (3.27) 

is used, where excip,£, = 0) = e^^^(p), ixcip,£, ~ 1) can be compute from 
the energy of the homogeneous electron gas in the configuration ^ = 1 and 
g(^) is an interpolation function, which is usually chosen as 

^^^^ 2 (2V3 _ 1) ■ (3-28) 
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3.4 Time— dependent density— functional theory 

DFT is a convenient theory to obtain the ground-state of a system through 
minimization of the energy as a functional of the density. However within 
DFT it is not possible to describe the evolution of a system perturbed by 
an external time dependent potential. The extension of DFT to the time 
domain, the Time Dependent DFT (TDDFT) relies on a formal extension of 
the HK theorem given by Runge and Gross (RG). In 1984 RG [21] showed 
that, given an initial time to where the system is in a state "^/q, is possible 
to establish a bijective correspondence between the time dependent external 
potential, the evolution of the density and the evolution of the wave-function. 
On the basis of this theorem they derived three schemes to calculate the time 
dependent density, one of which, a stationary action principle, can be seen as 
.he extension ot the second HKthc„,,,J3 . The theory provides a formalism 
where, at least in principle, it is possible to describe the evolution of a system 
even in the presence of strong perturbations. In this work we are interested in 
small perturbations and so we can restrict the analysis to the linear regime. 
The response function within TDDFT can be obtained through a single 
particle Schrodinger Equation, or Time-Dependent KS (TDKS) equation 
(that is the third scheme proposed in RG paper), by observing 

Spil) = XKs{l,l')5vKsa') (3.29) 

where xxsi^j'^) is the response function of the KS non interacting system, 
while SvKsi^) is the variation of the TDKS potential needed to follow the 
evolution of the density when an external potential (51^'^^* (1) is applied. From 
the previous section we know vks{^) = V^^^ll) + w/i[p](l) + Vxc[p]{l)- As- 
suming that the TDKS potential is identical to the potential of the static 



*The RG action functional has been found out to be wrong. A corrected functional 
has been provided some years later by Van Leeuwen 1221 . 
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theory. Its variation is, to linear order in the density, 



^...(1) = SV-^il) + '-^SpiV) + '-^^Spil'). (3.30) 

If we Compute the variation of the Hartree potential, defining the exchange 
correlation kernel fxc{^, 2) = Svxc[p]{i)/Sp{2) and remembering that Sp{l) = 
x(l, 1')6V'^^*{1') we obtain a Dyson equation for the response function: 

2) = XKs{l, 2) + XKsil, 1')/h.c(1', 2')x''(2', 2). (3.31) 

Here /i/a;c(l,2) is the sum of the Hartree plus the xc kernel. The most 
common approximations adopted for v^c and so for f^c are the adiabatic 
extensions of the ground-state approximations. However within TDDFT 
excited states, that can have different spin configurations, are described. 
For this reason it is convenient to derive the approximation from the SDFT 
energy functionals: 

r.c r , (5<^[p,C](xi,^i) , '5^g^-^[p«.,6J(xi) 

fa,aAp,^i^ltl,X2t2) = . , ^ - <5(il -^2) TTTT"^ ' 

0Pa2{^2,t2) dpa\[:>i2j 

(3.32) 

the variable ^ is eventually set to zero after the derivatives have been com- 
puted, for spin unpolarized systems to have a functionals of only the density. 
As for the ground-state, in the case of spin polarized systems TD-Spin-DFT 
(TDSDFT) can be used used, that is ^ is not set to zero. In particular start- 
ing from the LSDA we have the Adiabatic LDA/LSDA (ALDA/ALSDA) 
approximation 

Jxc (1,2) = (5(xi,X2)(5(ti,t2) 



dpa^ 

de'iriP.O , d^4r{p,0 \ (3 33) 
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Figure 3.2: Dealing witli tiie interaction. Form ref. [B]. 



Similarly to the case of DFT, TDDFT must be formally extended to TD- 
Current-DFT and/or TDSDFT if external magnetic fields are considered. 
However practical calculations indicate that this makes little difference in 
practice. 

3.5 The electron— electron interaction 

The electron-electron interaction is the main ingredient of Many Body sys- 
tems. Both DFT and MBPT offer a practical approach to deal with it. In 
DFT, and in particular in the KS scheme, the starting point is a system 
of not-interacting particles where the effect of the interaction is included 
through an external potential. In MBPT the classical part of the interac- 
tion, the Hartree potential only, is a local potential. Correlation is included 
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in the self-energy S which is an highly non local frequency dependent and, 
in the quasi particle picture, complex function. Compared to KS electrons 
QP can be seen as weakly interacting dressed electrons characterized by a 
finite lifetime. 

The main advantage of DFT compared to MBPT is that it's much less 
demanding from a computational point of view. This is the reason why nowa- 
days the majority of codes solve the MBPT equations starting from DFT 
calculations. Still the main drawback of DFT is that it does not offer a recipe 
to construct systematic approximations. For this reason attempts to improve 
the LDA are not always satisfactory. For example the Generalized Gradi- 
ent Approximation (GGA) does not always improve upon LDA in lattice 
constants and can even predict less accurate results [23l El] ■ Improvements 
to the LDA arc also derived from other approaches such as hybrid func- 
tional defined from the Hartree Fock (HF) approximation. Within these 
approaches DFT can include more physics at the price of loosing, at least in 
part, its main advantage. As a matter of fact the more physics is included 
beyond the LDA the more complicated are the functionals constructed. Fol- 
lowing this path, though many important results have been obtained, DFT 
is bound to become at a certain point more expensive than other competing 
approaches. 

In the present thesis we prefer to focus on MBPT in the first part where 
we will look for development of new approximations, in particular for the 
description of double excitations, and to look at LDA/DFT as the zero-order 
approach. For this same reason we will use DFT in the second part where 
we will explore the predictions of ab-initio calculations on carbon nanotubes 
in a magnetic field. 
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Part II 

Double Excitations 



Chapter 4 

Introduction to the 
problem 

To observe an object with our eyes we need to shine light on it. An absorp- 
tion experiment is nothing but a detailed study of what happen when the 
light hit the object. In this thesis we focus on the theoretical description of 
spectra in the visible energy range: the goal of theoretical spectroscopy is to 
give an accurate description of the microscopic process involved. The theory 
which describes the evolution of the many body electronic wave-function 
interacting with the electromagnetic field is the quantum electrodynamics. 
Unfortunately it is not possible to solve the exact equations and, as we have 
explained in the first part of the present work, approximations need to be 
introduced. Besides mathematical rigour, physical intuition can give essen- 
tial guidelines in the development of approximations. Relying on concepts 
emerging in the macroscopic world, obtained as the classical limit of quan- 
tum mechanics, it is often possible to give an intuitive interpretation of the 
microscopic world. As an example we can consider the idea of quasiparticles, 
which can capture many of the features observed in the absorption exper- 
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iments describing the main peaks as excitons. Similarly we have the idea 
of collective excitations, which is used to explain other possible features as 
plasmons. It is somehow surprising that, using these concepts, one gets an 
accurate description of the complicated microscopic time evolution described 
by the equations of quantum mechanics and quantum electrodynamics. 

However some of these simplified pictures fail in some situations. For 
example the description of the satellites discussed in Fig. p.3p requires a 
more elaborate representation then the single-particle picture. 

The existence of double excitations (DEs) is, in this sense, another situ- 
ation where an intuitive picture, like the one of quasiparticlc or the one of 
collective phenomena of the many body problem is problematic. It is still 
possible to describe the process as a result of the interaction among elec- 
trons: the light source hits an electron which is excited to an unoccupied 
state leaving a hole in the system. Due to the presence of both the hole and 
the extra electron the system can react by adjusting its configuration and ex- 
citing a second electron to an unoccupied state. This description is, however, 
not completely correct as it does not describe correctly the real evolution 
of the system. The alternative mathematical view is pictured in Fig. ()4.ip . 
The frequency-dependent kernel describes the correction to an independent 
particle picture due to correlation effects. The structure is very close to the 
one described in Fig. (|2.3p . The frequency dependence of the kernel reflects 
the fictitious time evolution used in MBPT to make the system evolve from 
a nn-interacting to an interacting eigenstate. DEs are then virtual processes 
needed to describe some features of absorption experiments. 

It is important to observe, however, that Fig. ()4.1|) describes a different 
process from the one of Fig. ()2.3p . In Fig. (j2.3p we are moving from an 
independent particle picture to a QP picture with a finite lifetime which, 
in the sense shown in Fig. (j2.4p . contains more physics than the simple 
representation in terms of real poles. The common concept, however, is the 
failure of the static approximation which does not give the correct number 
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Figure 4.1: The frequency-dependent kernel of the Dyson equation for the re- 
sponse function sphts the nn-interacting pole u)q in two poles u)a and ujb- The pole 
uj is obtained using a static kernel. From Ref. [25] . 

of solutions. This happens when the kernel of the Dyson equation has a pole 
at an energy close to an independent particle solution. 

The most common and widely used approach to the problem of DEs 
is given by post-HF methods based on the HF starting point which do 
not include correlatior0. However, while this is often reasonable for small 
molecules, in long ID molecular chains the effects of correlation are cru- 
cial |26| . On the other hand the currently used approximations to electron 
correlation in the state-of-the-art approaches to the description of optical 
excitations, TDDFT and MBPT, fail to capture the physics of DEs pTll^ . 
As post HF methods are designed to work in a small correlation regime, 

^ Correlation is due to the effects of the interaction in a quantum many-particles sys- 
tem. In the present thesis I use the following "language": if it's possible to save the HF 
interpretation as a first order approximation the correlation is considered small, if not the 
correlation is considered big. When DEs are important, the interpretation in terms of 
single particles (or HF orbitals) starts to break but not completely and indeed diagrams 
one order beyond the HF (or TDHF) are enough to include this effect. 
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their extension to realistic nano-structured materials is very demanding, 
if not practically impossible and solutions within the TDDFT and MBPT 
framework are needed. 

In the TDDFT approach, the excitation energies of a system are obtained 
from the nn-interacting Kohn-Sham (KS) cigenenergies solving Eq. (|3.3ip 
where the xc-effects are cast in the unknown xc-kernel /a;c[/5](r, r', i — t'). 
Most of the success of the scheme is due to the success of the ALDA, which, 
despite being extremely simple, is surprisingly accurate in the case of many 
isolated systems. Nevertheless the ALDA suffers from some deficiencies that 
cause TDDFT to fail in some cases, such as in the description of excited 
states with multiple-excitation character. The source of this failure has 
been traced back in the literature to the adiabatic approximation [271 US] , 
which neglects the frequency dependence of the true xc kernel: it turns out 
that it is precisely this frequency dependence of the kernel that takes into 
account all the many-electron excitation effects. 

In MBPT the neutral excitations of the system are obtained by solv- 
ing Eq. ([^75)1 . Similarly to TDDFT, xc effects in the BSE are cast in the 
four-point kernel S(l,2;3,4), which, unlike in TDDFT, can be written as 
a perturbative expansion. In the most common and widely used approxi- 
mation to this kernel, introduced in Eq. (|2.74[) . the xc effects in the BSE 
arc described by the screened Coulomb interaction W , which is considered 
static, thus ruling out the possibility of describing DEs. 

In the TDDFT literature several solutions to the double-excitation prob- 
lem have been proposed [551 HIl HH 1301 ISI] ■ Wang and Zeigler used a non 
collinear representation of the xc kernel [30] , which could be used to describe 
double-excitations, but only starting from the appropriate reference excited 
states, some of the weaknesses of this method are discussed in Ref. [5T| . 
Casida proposed a xc kernel which goes beyond the adiabatic approximation 
constructed from a superoperator formalism [2^ that contains as a special 
case the "dressed-TDDFT" recipe derived by Maitra et al. [53] . The dressed- 
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TDDFT approach however is not predictive since the existence of DEs must 
be defined "a priori" ; only very recently Huix-RotUant and Casida [32] pro- 
posed an extension of the dressed TDDFT method, which clarifies the rela- 
tion between Polarization Propagator (PP) approaches and the BSE meth- 
ods and is presently being tested on an extensive set of molecules |32] . 

The most commonly used state-of-the-art approaches suffer, hence, from 
different types of pathologies. On one side, post-HF methods, based on the 
uncorrelated HF scheme, are in general designed to describe isolated systems 
with the idea that the interactions among particles can be treated perturba- 
tively. In these, approximations are obtained by truncating the perturbative 
expansion to some finite order but always respecting key principles of quan- 
tum mechanics such as quantum statistics and Pauli exclusion principle. On 
the other side, in the BSE and TDDFT approaches correlation is treated to 
all orders of perturbation theory, but a well established method to include 
DEs does not exists yet. 

In the present chapter we provide a detailed description of the phenom- 
ena where DEs play a role and we illustrate which direction can be followed 
to include DEs in TDDFT and BSE. In particular we write the equations of 
both approaches in a common formalism (see Eq. (|4.4p ) and we show that 
DEs are not described in the standard approximations because, in both, the 
kernel of Eq. (|4.4p is taken static. We focus on the BSE scheme, because 
within the MBPT it is more straightforward to look for practical approxima- 
tions, in order to introduce a frequency-dependent kernel. First we solve the 
mathematical problems which arise trying to construct the dynamical ver- 
sion, i.e. with a dynamical kernel, of Eq. (|4.4|) and then we show that indeed 
the frequency-dependent kernel obtained relaxing the static approximation 
in the standard BSE scheme capture, at least in part, the effects of DEs. The 
inclusion of this result within TDDFT is obtained thanks to the common 
language established, following an alternative, but similar, procedure to the 
one of Ref . [33] . 
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Only in the following chapter, we will finally construct a fully consistent 
approximation to the problem of DEs within the BSE scheme. The explicit 
connection with TDDFT is not explored anymore in Ch. [S] the extension 
of the results obtained to TDDFT, following the method illustrated in the 
present chapter, could be a possible development. 



4.1 Double excitations in quantum chemistry 

The concept of DEs is well known in the quantum chemistry literature. 
The absorption spectrum of a system is described with methods that, in 
contrast to DFT and MBPT, arc based on the many body wave-function. 
The starting point is usually the Hartrce-Fock (HF) approach, where the 
many body ground state is approximated with a single Slater determinant 
of one particles wave-functions: 



HF 



'0l(xi) ■(/'l(x2) ■•• -iAlCXn) 
V'2(xi) ■(/'2(X2) ... '02(Xn) 

V'„(xi) V'«(X2) ... -iAnCXn) 



(4.1) 



The excited state wave-functions can be expressed as a linear combination 
of Slater determinants, each of which is related to the HF ground-state 
through some excitation operator. For example, considering only single- 
particle excitatins one can write: 



HF 



(4.2) 



which we refer here as with / labeling a particular excited-state. In 

the language of MBPT the HF approximation is equivalent to approximating 
the self-energy to its first order in the bare interaction. The approximation 
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Rh-h [^-^-l R^.^ [a.u.] 

(a) states with Su symmetry (b) states with Sg symmetry 

Figure 4.2: The lowest energies excitations energies of the H2 molecule as a 
function of the distance among the nuclei. Solid lines are CI results, while dotted 
and dashed lines are TDDFT calculation with two different functionals which does 
not include DEs. The excited state ^E^ mainly composed by a DE is marked with 
dots. From Ref. [34l- 



introduced by Eq. (j4.2p . instead, is equivalent, in the linear response regime, 
to approximating the BSE kernel to first order in the bare interaction within 
the Tamm-Dancoff Approximation (TDA) which neglects the coupling of the 
excitation space with the de-excitations one. The scheme obtained relaxing 
the TDA is known as Time Dependent HF (TDHF). 

It is then natural to relax the introduced approximations to include single 
and higher order excitations in both the ground state and the excited state 
wave-function. This is the so called configuration interaction (CI) expansion: 

ij ijhk 
ijhk 

In practice the full CI approach can be hardly used for a number of elec- 
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trons larger than ICo, because the exponents grown at the number terms 
needed. However CI offers a clear mathematical description of multiple ex- 
citations, including DEs. 

In the next chapter we will present an alternative Post-HF method, called 
second Random Phase Approximation (sRPA) which we will use to investi- 
gate more in detail the connection between the MBPT and the DEs concept. 
Here we want to start by using the CI scheme to understand how DEs are 
identified. We will also present two examples where the DEs have been iden- 
tified to play an important role using CI calculations. The first example is 
constituted by small molecules, like the H2 molecule considered here, for 
which DEs are known to be important in the description of molecular disso- 
ciation; the second is constituted by polyenic carbon chains saturated with 
hydrogen atoms. 

H2 is a very simple molecule and DEs play a crucial role already for such 
a small system. In Fig. (|4.2p we can see that during the dissociation pro- 
cess even the qualitative behavior of the excitation energies is wrong within 
TDDFT, if compared with a virtually exact full CI approach. In particular 
for the states with g symmetry, panel &, we see that within TDDFT we have 
one excited state less than within CI. The double excited CI solution start at 
high energy at the equilibrium distance but then decrease in energy during 
the dissociation process becoming the lowest energy solution. 

For polyenes, on the other hand, the description of DEs is known to be 
important already at the equilibrium geometry. The structure of the chains 
is represented in Fig. (|4.3p . together with a plot of the contribution of DEs 
to the ground state and the first three excited states energies as a function 
of the chain length. The contribution to the ground state increases with the 
dimension of the chain from less then 10% to around 25% while for two of the 
three excited states is above 50% (reaching the 75%) even for the shortest 



^CI in the "singles and doubles" approximation (CISD) has recentrly been performed 
for systems up to almost 100 electrons, though under some approximations 135 1 . 
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(a) polyenes structure (b) Double excitations 



Figure 4.3: Structure of a polyene chain, (a). The fundamental element of a 
chain is represented between brackets. At each vertex a carbon atom saturated 
with hydorgen is present. The weight of doubly excited determinants in the wave- 
functions of the electronic ground state 1^A~ and the energetically lowest three 
excited states i^A~, 1^A~ and 1^A~ are given in percent. From Ref. [36] . 



chains. 

Theoretically, DEs arc predicted to be important for any open-shell 
molecule and in general for any system where the energy of a DE can be 
degenerate, or almost degenerate, with the energy of a single excitation, as 
shown in Fig. ()4.1|) . For open shell systems, in particular, the inclusion of 
DEs is imposed by spin symmetry requirements which forces single excited 
configuration to mix with double excited ones. This is shown in Fig. (|4.4[) for 
a very simple three electrons system. The four single excited configuration 
represented must be mixed among themselves to obtain an eigenstate of the 
total spin operator 5'^. However such operator contains a term of the form 
al^aj^Sj-t-dij, which flips the spin of two particles. Applying these terms to 
the configurations ip^ = al^ai■f^po or ■04 = aj^^a^^o a doubly excited config- 
uration -05 = a^^Oj^aj^diit/jQ is obtained. This is needed to construct the 
excited state wave-function with the correct spin symmetry |37[ 127] . 
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Figure 4.4: Excited states in a spin polarized model with three electrons and 
three levels. The ground state configuration is plotted on the top, while the doubly 
excited configuration is highlighted with respect to the others. From Ref. [27] 



Satellites in the absorption spectrum 

We have already drawn a connection between the concept of satellites in the 
photoemission spectrum, illustrated in Ch. [2] and the mathematical defini- 
tion of DEs. In fact DEs can be seen as satellites in the absorption spectra: 
satellites in photoemission spectra are generated by the frequency depen- 
dence of the self-energy, whereas DEs by that of the BSE kernel. However, 
as opposed to photoemission spectra, the concept of satellites (due to pure 
many body effects) is not used for absorption spectra in the visible/UV 
range. The main reason being probably related to the different structure of 
the self-energy and of of the BSE kernel. The Self Energy is in fact usu- 
ally almost diagonal and as a result the independent particle peaks are only 
shifted without involving a mixing of the independent particle transitions. 
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(a) Theory: TDDFT and LDA+U (b) Experiment: NIXS 



Figure 4.5: Non- Resonant Inelastic Xray Spectrum (NIXS) of NiO solving a series 
of tightly bound Frenkel excitons. The experimental spectrum [38] shows more 
peaks than theoretically predicted [35]. In the work on Lee et al. [35] a frequency- 
dependent kernel is proposed as a possible solution to describe the experimental 
peaks (indicated with three arrows in (a)). Within the adiabatic description only 
two peaks are obtained. 

On the other hand the kernel of the BSE is strongly non-diagonal and an 
exciton is composed by different independent particle transitions. As a con- 
sequence, usualljlfl, in the QP spectrum it is easy to isolate a peak which 
is not related to any independent-particle (IP) transition and so to trace it 
back to a satellite. On the other hand in the absorption spectrum many IP 
transitions usually mix in different ways to give well-defined cxcitonic peaks 
together with other less intense and dark peaks. Many features are already 
present in the spectrum. The inclusion of DEs can shift the excitonic peaks 
and increase the richness of the spectrum, however it is difficult to isolate a 
DE satellite without a Cl-like analysis of the peak composition. 

NiO has been shown to be a good candidate system where DEs can be 
found in the absorption spectrum. In two recent works the Non-Resonant 

•^In strongly correlated materials the QP peak can be difficult to isolate as well. 
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Inelastic Xray Spectrum (NIXS) of NiO has been investigated both theoreti- 
cahy and experimentally [511 ISH] ■ At low energy NiO presents tightly bound 
Frenkel excitons (see Fig. (|4.5|) '). However the theoretical description within 
adiabatic TDDFT is not able to reproduce the detailed structure of the spec- 
trum (the theory reproduces only 2 of the 3 peaks observed experimentally). 
In their work, Lee et al. [32] indicate the need to go beyond the adiabatic 
approximation for a correct description of the spectrum: «...to allow any 
fine ("multi-plets") structure in strongly interacting systems, a non- adiabatic 
kernel is absolutely necessary (■■■). Obviously, this is one key aspect that al- 
most all the existing approximate functionals lack and presents an essential 
and necessary step toward a proper description of local excitations in strongly 
interacting systems, within all the existing theoretical frameworks. >>0. 



Double excitations, in which basis set? 

Previously in this chapter we have defined DEs from the CI expansion in 
terms of the HF wave-functions. However some of the experimental evidence 
we have illustrated has been compared against TDDFT calculations. It is 
then legitimate to question if we can use the same definition of DEs using 
the KS wave-functions, which, unlike the HF ones, do not have necessarily 
any physical mcaninjfl. 

As we stated in the introduction, when a feature of the spectrum need 
to be described in terms of DEs, the single particle description is breaks 
down. Then the HF wave-functions lose their physical meaning and have to 
be regarded only as a basis-set for the expansion of the Many-Body wave- 
function, and the same holds for KS wave-functions. The open question is 

■^However alternative explanations to the extra peak observed in NiO, whieh do not 
make use of the concept of DEs, already exist in the literature 1401 . 

^There are several articles indicating that KS orbitals are typically good approxima- 
tions to Dyson orbitals 1411 1421 |43I . However this is not true for any system. In corre- 
lated materials the overlap between KS wave— function and QP wave— function has been 
proven 
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if a double excitation in a basis set can be described as a single excitation 
in another basis set. In fact the single particle transitions, either KS or HF, 
do not form a complete basis-set in the space of the many-body excitation 
operators. Indeed the two span two different spaces. 

As a partial answer it is possible to consider the number of excited states: 
the inclusion of DEs gives more solutions than the number of single particle 
excitations initially considered. However the space of single particle transi- 
tions is infinite and so it is the number of solutions. 



4.2 Double excitations in the many body ap- 
proach 

Both TDDFT and the BSE provides exact equations for the description 
of absorption spectra, therefore they can, in principle, describe DEs. In 
practice, however, approximations to the many-body effects of the system 
are needed, and the currently used ones, namely ALDA for TDDFT and the 
statically screened interaction within BSE, fail to reproduce DEs. This is 
why in this work we go beyond the standard approximations and we derive 
a new kernel for BSE and TDDFT that can properly take into account DEs. 

In the following we will introduce the single-particle transition space and 
we will show why a static kernel cannot describe DEs. By relaxing the static 
approximation to the kernel one can get DEs. We will illustrate this using 
the BSE, where the kernel has a clear physical meaning and includes, in a 
natural way, the many-electron excitations of the system. From the BSE 
kernel one can then obtain the TDDFT kernel using the technique of Rcf. 
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The Dyson equations in the transition space 

The main advantage of writing the equations for the response function in 
the space of single-particle transitions is that it offers a clear interpretation 
in terms of the single-particle wave-function. In the present work wc will 
assume that the differences among the QP and the KS wave-functions are 
small and we will project the TDDFT and BS equations in the KS basis set. 
This assumption is usually done in the implementation of the BSE scheme 
within many ab-initio codes and it is justified in many materials though 
exceptions have been found [TT] . In particular the starting point is usually a 
DFT-LDA calculation with self-energy effects introduced according to the 
scheme outlined in Ch. [21 Then the excitonic spectra can be computed 
within TDDFT starting from the V^f and ef^"^ or from the MBPT-GW^ 
approach where the BSE is solved starting from i/jf^ and ef^. 

We start from a Dyson equation for a generalized four-points response 
function in the frequency domain 



i(l,2;3,4|w) = Lo(l,2;3,4|tj)+ 

Lo(l, 2; 1', 2»A-(1', 2'; 3', 4'|0)Z(3', 4'; 3, 4\lu), (4.4) 



where we have already introduced the static approximation for the kernel 
K. Eq. (|4.4p can be obtained from Eq. (|2.73|) by using 




Lo(l,2;3,4|cj) 



/ 



d{t2 



^*^-*i)Lo(lti,2ii;3i2,4i2), 




S{l,3)S{2,4)W{lti,2t2)S{t2 - t 
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we have introduced L to be distinguished from L. L is a function of four 
space variables as L but two time variables only. An equation of the form 
of Eq. (|4.4p can be obtained from Eq. p.3ip using 

i(l,2;3,4|c^) = j dL.e-'"(*^-*i)x(lii; 3*2)5(1,2)5(3, 4), 
io(l, 2; 3, 4|w) = j dwe-*"(*^-*i)xKs(l*i; 3*2)5(1, 2)5(3, 4), 



i^(l,2;3,4|0) j duje-"^^*'-*'^ (^w{l,3) + 



5pti(3) 



5(1,2)5(3,4). 



Notice that not only the kernel is different, but also the four-points functions 
L and Lq differ in the two cases. In the space of single particle wave-functions 
the equation fori reads! 

Uj^uj) ^ ,,fc(a;) + L%^,,j,{uj)K,,j,^H'k-{Q)Lh'k'M^)^ (4-5) 
where we used the following change of basis 

Uj,hk{^) = / rf^xi...d^X4 ^»(xi)V'*(x2)L(l,2;3,4|w)V'h(xi)^/.fc(xi) 

J (4.6) 

= Z(1,2;3,4|l.) \hk). 

Here i is a generalized index for the KS wave-function containing all quan- 
tum numbers, spin included. Within this basis-set the generalized response 
function Lq is diagonal 



?o _ 5t,fc5j,/,(/j - fi ) 
w - (e* - ej) 



L%,Hk = '^' 'T' (4.7) 



^In the TDDFT formalism, in order to obtain simpler equations, when moving to the 
single-particle wave— functions space, the delta functions are discarded in the definition of 
Lq which is constructed from Lxs(l, 2; 3, 4) = G(l, 3)G(4, 2). Only in this way Lq. 
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We have dropped the ir] factors in order to avoid having two different expres- 
sions for the retarded and the time-ordered response function. If needed one 
can restore them remembering that we deal with retarded quantities within 



TDDFT and with time-ordered quantities within MBPtH Moreover we 
have introduced the occupation factorial, = 1 if i is occupied and fi = 
otherwise, to have a compact expression for Lq. Using Eq. (|4.7p and writing 

^0 

L^^ solving an eigenvalue problem: 



the matrix equation in the form L ^ — L^^ ~ K one can find the zeros of 



H^PAi = LUiAi, (4.8) 



with HijM = {ei - ej)Si^kSj,h + \fJj~~Ti Kij^hk Vfh - Ik- Taking explicitly 
into account the occupation factors the eigenvalue problem can be recast in 
the electron-hole {eh) and hole-electron [he) transitions (see [H] for TDDFT 
and [ig for BSE): 



coup ^res 





= ^i{ ^ . ]{ J ]■ (4-9) 



We defined the resonant term Hres = Heh.e'h' and the coupling term Hcoup = 
Heh,h'e' and we have assumed that A'(l, 2; 3, 4|a;) is real. Let us call 5* the 
matrix defining the metric of the system: 




(4.10) 

The block form of Eq. ()4.9p allows us to have a clear interpretation of the 



^One have to take care of this difference if a connection among the two theories need to 
be estabhshed. The difficulties which can arise trying to combine theories with different 
time ordering can be overcome thanks to the Keldish contour techniques. 

*The occupation factors can be used as Fermi distribution functions to introduce a 
numerical smearing to get a faster convergence. The same Fermi functions are sometimes 
used to introduce a temperature dependence in the equations. 
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physics involved. In particular the resonant part of the Hamiltonian de- 
scribes the neutral excitations of the system. In terms of a linear composition 
of single eh transitions. This can be directly compared with 'i'f^ defined in 
Eq. (|4.2|) . The term H*^^ is the anti-resonant part of the Hamiltonian and 
describes the de-excitations. Finally H*^^^ describes the coupling among 
the eh and the he space. When these are different from zero the many-body 
excitations contains terms involving single-particle de-excitations. These 
latter processes are clearly forbidden if one consider the HF ground state 
and for this reason the coupling terms are said to describe ground state 
correlation [15] , 

While the eigenvalues of the problem gives the excitations energies the 
eigenvector can be used to construct excitation operators. Since we are 
working within the DFT basis set, it is tempting to describe the excited 
state wave-function directly applying the excitations operator to the ground 
state, in analogy to HF-based methids. However one has to keep in mind 
that the TDDFT and MBPT linear response equations have been derived 
starting from a variation of the density with no assumption on the wave- 
function. For this reason the excitation operator should be applied not to the 
DFT wave-function but to the correlated Many-Body wave-function. Only 
the TDA, which assumes that de-excitations do not need to be considered, 
is consistent with the approximation of the ground state as a single slater 
determinant. Within this approximation excited states wave-functions can 
be constructed from the KS ground state and interpreted, as it is often done 
for KS the ground-state, as approximations to the real many body wave- 
function. However only within the TDHF scheme and starting from the HF 
ground-state the approximation is formally correct. 

Moreover, with the exact kernel, Eq. (|4.9p would give the exact excita- 
tions energies, but not the full excitation operators, since only the single- 
particle part is accessible by construction. This means that multi-particle 
transitions must be hidden in the kernel of the equation. It becomes now 
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clear why the eigenvalue problem ()4.9p cannot give DEs, unless the kernel 
is frequency-dependent. In the next chapter we will clarify these points 
exploring the Second-RPA method where an equivalent eigenvalue problem 
will be derived exactly from the projection on the space of single particle 
transition of an exact equation for the excitation operators. 



Eq. ()4.9p is equivalent to a Dyson equation for the generalized response 
function. From its solution wc can construct the response function x(l, 2|aj) = 
L{1, 1; 2, 2\lj). As explained in Ch. [T] we can use x to define the dielectric 
function {e^^{Lu)) = (1 +wx) and the polarizability (a) = (x^x). As an ex- 
ample we write here the esplicit connection of the macroscopic measurable 
quantity a(w) with the Eq. (|4.9p . 



The expression for the generalized response function in terms of the eigen- 
values and the eigenvectors of Eq. ()4.9|) is 




(4.11) 



Using Eq. (|4.1ip the polarizability can be written as 




ij.hk 




(4.12) 



The direct expression for the dielectric function in the space of transitions, 
instead, can be used to describe the concept of local fields effect if the basis- 
set of the block wave-functions is used. This is done in App. [X] 
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4.3 The dynamical Bethe— Salpeter equation 
(step I) 

In order to explore the effects of a dynamical kernel on the description of DEs 
as start we simply relax the static approximation to the screened interaction 
in the BSE §|. However wc immediately realize that, by using a frequency- 
dependent interaction W{lu), Eq. (|4.5p cannot be written any more as a 
simple matrix equation, as it involves a convolution in the frequency space. 
To analyze the problem we rewrite the exact BSE for L(uj) obtained by using 
a frequency-dependent kernel: 

L,,^M = L0.,,(^) + dLo'dLo"Ll^^,^,{u,u') 

^i'j',h'k'{^,^',(^")Lh'k'.hk{i^,i^"), (4.13) 

The kernel S depends on four time variables making it impossible to contract 
the variables in the last term of the r.h.s. of Eq. ()4.13|) for L and L^, as 
shown in Fig. (|4.6p in order to obtain a closed equation for L. This problem 
does not appear within TDDFT where we have only two point (and so two 
times) quantities. For the Fourier transform we have adopted the following 
conventions [47 1 l33j: 

L(w, uj', Lo") ^ J dTdT'dT"e'^^e'^'^' e'^"^" Liti,t[,t2,t'2) (4.14) 

with T = (ti + t[)/2 - {t2 + t'2)/2, t' = ti - t[ and r' = ta - ^2- 

Commonly the BSE kernel follows from the GW approach to the self- 
energy and it reads 

Sij,/ifc(a;" - uj') = Wijjih - Wij.hki^^" - w'), (4.15) 



^We choose BSE as starting point because within TDDFT wc do not have a straight- 
forward way to insert a frequency-dependency in the fxc kernel 
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Figure 4.6: The BSE with a static kernel (upper diagrams) is compared with the 
BSE with an exact kernel (lower diagrams). The contraction of time variables to 
obtain i — > L is represented with small dashed lines. The static approximation 
for the kernel is represented as a collapse of the time dependency to a single point 
highlighted in the figure. When using the exact kernel there is no time contraction 
and it is not possible to do the contraction L — >■ L on the last term of the r.h.s. . 
The BSE does not reduce to a closed equation for L. 

where we considered the full frequency dependence in the screened interac- 
tion. To obtain a Dyson equation we insert the identities LqL^^ and LL~^ 
on the left and of the right of the second term on the l.h.s. of Eq. (|4.13p and 
define a new kernel 

LO.,^-(o.,wOW.,-,.fe(u;'' -^')i°fc,hJ'^,c.'')LL,M. '\u;). (4.16) 

Here we have approximated L c^i Lq (linearization) in the kernel expression 
and we have used the fact that Lq is diagonal in configuration space. In this 
way we have a closed equation, the dynamical BSE (DBSE), that wc can 
project in the transition space in order to obtain an eigenvalue problem: 

U,,hk{^) = Llhui^) + L%,,.y{uj)k,,y ,h'k'{^)Lh'k' ,hk{^); (4.17) 

which is formally identical to SBSE except for the presence of a frequency- 
dependent kernel K{ijj) = S2(w). Note that the if a static kernel is used the 
DBSE reduces exactly to the usual Static BSE (SBSE). 



4.3 The dynamical Bethe Salpeter equation (step I) 



77 



From Eq. ()4.17p we can construct an eigenvalue equation with frequency- 
dependent excitonic Haniiltonian 

H^P{uj)Ai{u) =ui{u)Ai{uj). (4.18) 

Note that Eq. (|4.18p is similar to the frequency-dependent eigenvalue equa- 
tion obtained by Strinati in Ref. [4 7) . 



A first analysis of the kernel 

The frequency integrals in the definition of K{uj) can be performed analyti- 
cally as we know the frequency dependence of both Lq{uj) and W{ll!). As for 
the static screening, we have here the problem that the kernel of the equation 
should, in principle, depend on the solution of the equation itself. However, 
in the SBSE the exact position of the poles in the construction of the kernel 
is not important and LDA eigenvalues are used to construct the screening 
instead of the exact poles of the response function or of the QP energies. 
QP eigenenergies are not chosen in order to prevent the larger QP gap to 
underestimate the screeninJ"!. Indeed LDA gap is often comparable with 
the optical gap due to partial cancellation of self-energy and kernel effects; 
see Fig. (j4.7p . However we will show in Ch. ([5]) that the DBSE solutions 
are very sensitive to the exact position of the poles (see also Fig. ()4.1|) ). 
Therefore a KS screening might not be accurate anymore. To overcome this 
problem one can use RPA or static BSE energies. From now on we will use 
the RPA energies il^, with eigenvectors , which are solutions of Eq. (|4.9p 
with /v(l,2,3,4|0) = (5(1, 2)5(3, 4)w(l, 3). With this choice we can write 



^"in a recent work 1481 the use of QP energies in the construction of the kernel has been 
proven to give better results for some materials. Other work on this point can be found 
in [3911501 
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Independent particles 



GW 



BSE 



Figure 4.7: Self-energy corrections open the LDA gap due to a better description 
of the (screened) exchange effect which tends to push particles away. However the 
(screened) attraction between the electron and the hole gives a binding energy 
which partially compensates the self-energy effects and so the value of the optical 
gap is often closer to the LDA gap rather than to the QP one. From Ref. [51] 



the frequency-dependent kernel within TDA as 



where the indexes (nq) run over the possible eh couples only and where the 
factor 2 is obtained using R'^^q-^ = R^qn) which holds within TDA. 

A frequency-dependent kernel imposes to solve Eq. (|4.9|) self-consistcntly 
with respect to the frequency. In practice we need the solution of the equa- 
tion Ej{uj) = oj with Ej{u!) the eigenenergies of a DBSE at a given frequency. 
This assumption resembles the QP concept as it relies on the assumption 





''V(n'q')(ih)i 



(4.19) 



4.3 The dynamical Bethe Salpeter equation (step I) 



79 



that 



Ai^{u:)SijAj^{u:) 



Ai^{Ej)Sj,jA^,{Ej) 
(w - Ei) 



(4.20) 



For further details on this idea see App. [B] 

Connection with time dependent density functional the- 
ory 

In sec. 14.21 Eq. (|4.5p has established a common language between the SBSE 
and TDDFT. The DBSE, Eq. dHH]), is also the generahzation of Eq. (|T5|) 
to the frequency-dependent case. Indeed, with the choice K{lo) = fnxci^) 
in Eq. (|iT7)) . the TDDFT formalism is recovered. 

A kernel for the TDDFT can be derived taking advantage of this common 
language. Under the assumption that the differences between the QP and 
the KS wave-functions can be neglected we obtain 



where the term Xa's(^) ^ Xo^{<^) takes into account for the difference be- 
tween the QP and the KS starting point. The term S2(w) is the space 
contracted version of S2(w), obtained changing the four points Lq with the 
two points xo in the definition of the latter. Indeed, choosing the static ap- 
proximation for the screened interaction in the definition of ^2(0;), we obtain 
the so called Nanoquanta kernel [TT] . 

From now on wc will work in the BSE framework, focusing our attention 
on the problem of DEs. However, under the assumptions of this section, a 
TDDFT kernel can be derived using Eq. (|4.2ip and using a modified DBSE 



(4.21) 
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kernel (here ^2(0-') ) which can be derived as long as the spatial contraction 
Lq — > Xo is possible. This will be the case for the kernel proposed in the 
next chapteio. For a further study of the problem of the construction of 
a TDDFT kernel from a many-body formalism we address the reader to 
the references [3U [H] . In particular the spatial contraction (or localization) 
is extensively discussed in [35]. As pointed out there however, it is not 
necessary if the only goal is to obtain excitations energy, since one can do a 
DPT based BSE calculation. This simply because x ^^'^ L have the same 
poles; on the other hand if one is interested in intrinsic TDDFT quantities, 
as for example the time evolution of the density, the localization process can 
be crucial [35]. 

Some preliminary tests on a model system 

The performance of the DBSE can be tested in a two electrons and two levels 
model, that is the simplest possible system where a DE can appear. We will 
work here within the Tamm-Dancoff approximation to keep the equations 
as simple as possible. We look for the solution of the eigenvalue problem 
solving the equation det{H^P{oj) — Icj) = 0. For our model we obtain: 




(4.22) 



where we defined V = Wyc,vc , W{i^) = (^2)^,^ „^ {^) , Ae = Cc — e« and we 
used the fact that W{uj) is diagonal in the spin space. The RPA solutions 
Qi^2 needed to construct the W{uj) term have eigenvectors i?i = l/\/2 (1 1)* 
and R2 = I/V2 (1 - 1)*, thus we get 



W{uj) = A + 



B 



(4.23) 



cj — il 



1 



However to define the contracted kernel the substitution Ls — > Xs instead of Lq — ^ XO 
will have to be considered. Ls will be defined in Ch. [5] 
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where with A = Wy^^cc and B = 2Re[wccvcWvcvv]- Eq. ()4.22p has four solu- 
tions, altliougli one expects only three for this system, i.e., a singlet single 
excitation, a triplet single excitation, and a singlet DE. One of the four 
is an unphysical state. We argue that the occurrence of this extra pole is 
related to the self-screening interaction that the GW approximation to the 
self-energy suffers from. This is related to the fact that W is the test charge- 
test charge screening, whereas the charges to be screened arc fermions, not 
classical charges. This can be cured by introducing a vertex correction to the 
self-energy. Indeed, if one considers only one electron in this model system, 
then Eq. (|4.18p produces two poles, one corresponding to a single excitation 
and the other one, unphysical, corresponding to a DE. In this case, there 
are no dynamical self-energy effects involved, and the extra pole arises, in- 
deed, from the fact that the electron screens itself. Wc can recognize the 
spurious solution by solving Eq. (j4.22p independently of the dynamical struc- 
ture of W. We then obtain two groups of solutions: one for singlet states, 
oj = Ae + 2V — W, and one for triplet states, w = Ae + 2V — W. Since the 
excited state involving a DE is a singlet, the correct double-excitation energy 
is the one coming from the singlet-group solutions. The four solutions {(^1^2 
the singlet solutions, and a;3.4 the triplet solutions) are 

2Ae + fii - A T J{ni-2V + A)^ - iB 

UJl,2 = 7) 

4 . (4.24) 

2Ae + ni- At + A)'^ - AB 
U^3A = ^ 

The energy W4 (the solution with the sign +) is a spurious pole. In the next 
chapter we will understand better the origin of spurious excitation energies 
and we will show how to derive an approximation to the DBSE kernel which 
does not suffer of this problem. 
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Chapter 5 

A new approach 



As we have shown in Ch. 21 DEs arc essential for the description of the 
optically excited states in open-shell molecules |27j : however they can play 
an important role also in closed-shell systems, such as in polyenes, where 
the lowest-lying sing let state is known to have a HOMO^-LUMO^ double- 
excitation character [28] . The theoretical description of Double Excitations 
(DE) in conjugated polymers constitutes an important challenge for the 
state-of-the-art approaches used in physics and physical chemistry. 

On one side there are the Post-HF methods that descibes DEs in a 
natural way, but at the price of a very demanding description of correlation 
effects; on the other side there are methods as TDDFT and BSE, which 
treat better the correlation, but within the standard approxmiations cannot 
capture the physics of DEs. The limitation of the latter approaches lies in 
the adiabatic approximation to the exchange-correlation effects. 

In Ch. 21 we showed that simply relaxing this approximatin, DEs are in 
fact described; however, together with the desired excitations, non-physical 
excitations also appear. Spurious excitations have been interpreted as due 
to the self-screening error embodied in the GW self-energy [331 ■ 
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In this chapter we investigate more in details this problem showing that 
uncontrolled effects, such as unphysical excitations, can appear as the quan- 
tum statistics and Pauli exclusion principle are easely broken in simple ap- 
proximations, like the one we introduced. For this reason here we propose 
a novel approach to describe DEs in correlated materials by embodying the 
mathematical properties of post-HF methods (here we will use as reference 
the second-RPA) in a coherent Many-Body framework. In order to achieve 
this we first define the conditions for a Number Conserving (NC) approach, 
which avoids the appearance of spurious excitations; we then embody the 
NC condition in an extension to the BSE that describes DEs in a consistent 
manner. 



5.1 The second random phase approximation 

A number— conserving approach 

The second-RPA (sRPA) is a particular appealing starting point because the 
scheme is directly derived approximating the many body excitation operators 
to DEs [531 : 




where dj / d^ are creator / annichilation operators in a single particle wave- 
function basis set, the HF wave-functions are used in the original derivation 
of the sRPA equations, and w/ are the excitation energies. The scheme can 
be constructed inserting Eq. (j5.1|) in a double commutator equation which 





(5.1) 
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is satisfied by the operator Oi, [S3]: 

{HF\ \r, \H,dl\] \HF)=u, {HF\ [i?, 0^1 \HF). 



(5.2) 



Here \HF) is the HF ground state, H is the many body Hamiltonian and R 
is an operator in the same space of the excitation operators 0,y. The result 
can be written in the form of an eigenvalue equation: 



A 

-B* 



where 



A 



and 



A 



ij,hk 



A, 



A. 
A,, 



B 

-A* 



B 



(5.3) 



A. 



ij,hk 



■jmn^hk -^ijmn^hkpq 



A. 



ij.hkpq 



(5.4) 



Mjmnjik -^ijmn^hkpq 



X 



(1) 



(1) 



ijmn / \ tjnin / 

The elements of A are obtained from Eq. (|5.2p using Eq. ()5.ip . for example 



A. 



ij,hkpq 



(HF\ 



\HF) 



(5.6) 



and similarly for the other components. The elements of B have a similar 
form as the elements of A, the only differences being: (i) in the operators 
on the right of the Hamiltonian particle-hole indexes are inverted, (ii) there 
is a minus sign. 

Eq. ()5.3p has the same formal properties of the RPA-TDHF equations and 
this guarantees that spectral sum-rules are respectecf^. Moreover measurable 
quantities can be constructed from the solution of the problem using the 



^ While the spectral sum rules are respected, the excitation energies are well described 
but, often, the oscillator strengths are not 1321 . 
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same equation of the RPA method [55] , 

Second random phase approximation and the folding 

The Hamihonian associated to the sRPA equation of motion can then be 
written [331 \SM ^he Fock space of single and DEs 



Here S and D represent, respectively, the Hamiltonian in the space of sin- 
gle excitations (dimension A^^ x Ng) and of DEs (dimension Nd x N^)- C 
represents the coupling between single and DEs. The number of eigenvalues 
of Eq. ()5.7p is, thus, Ng + Nd- ei and 62 are the sRPA excitation operator 
components [S31 [SI] in the singles and doubles subspaces, respectively. 

The question now is how to obtain these poles working only in the 
space of single excitations, without introducing explicitly the doubles sub- 
space. This step is crucial to create a link between the sRPA, Eq. (|4.12p . 
and the BSE, which is strictly defined only in the singles subspace. To create 
this link we fold the total Hamiltonian matrix in the R^" subspace [S31 [SI] . 
This is done by expressing 62 in terms of ei, and then solving the equation 
for ei: 



with E{uj) = C{uji - Dy^C''. Eq. ([ST]) and Eq. ([EH]), then, have the 
same Ng + Nd eigenvalues but Eq. ()5.8p is solved in the single-excitation 
subspace, and the frequency-dependent kernel S(a;) takes into account the 
down-folding of the double-excitation space to the single-excitation space. 
The correct structure of the S kernel is thus crucial to get the correct number 
of solutions. In particular, if D can be diagonalized, then Eq. (|5.8p can be 




(5.7) 



{S -|-^(a;))ei ^w/ei, 



(5.8) 
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written in terms of ttie diagonal matrix D' — DU: 



witlr K = C'C't and C" = CU. 

The explicit expression for the S kernel of sRPA can be obtained [53l [54] , 
within the TDA, starting from Eq. (|5.2p and Eq. ()5.ip and constructing all 
the matrix elements A and B: 

„ , . ^ii0)Anm){pq)Cl^rn)(pq),{hk) , . 

^(^JUhk){^)= -— 77-3-— T-^TT' ^^-10) 



, w , w - (e„ - e„ + - eg) 

(ng)(mp) 



with 



- {ri ^ p} - {m O g} + {(nm) O {pq)]^ ■ (5-11) 

Here, in the space of DEs, the matrix elements of the interaction term in the 
hamiltonin has been neglected; are the poles of the HP one particle GF, 
Ghf, whereas 

'"(i]),{hk) ^ J dxdx'(j)*{x)(l),{x)v{xx')(l)k{x')(j)l{x'), (5.12) 

are the projections of the Coulomb interaction in the space of single-particle 
wave-functions. The structure of Eq. ()5.10p is the same of the kernel in 
Eq. (|5.9p . A key property of the S kernel is that it is unchanged under — {n -H- 
p}, —{m <-> q} (Pauli exclusion principle) and {{nm) <-> (pq)} (particle 
indistinguishably) transformations due to the symmetry of the Ci^ij)^(^nm)(pq) 
factors. 
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Therefore the algebraic structure of Eq. i5.10\) ensures the respect of the 
particle indistinguishability and of the Pauli exclusion principle which con- 
stitute necessary conditions for a number-conserving (NC) theory of DEs. 

This can be shown in detail by solving the characteristic equation of the 
eigenvalue problem Eq. (|5.8p . i.e. detiuj — 5 — = 0. Using the non 

linearity of the determinant operator, 

, , , det(K^i'>) 

deti — = ^— 5.13 

where Ng is the dimension of the matrix K, and ^ stands for the set of 
indexes {(nm){pq)}, and exploiting the relation (55) 

det{A + B) ^ minor (A) minor (B), (5.14) 

Pr,Pc 

with Pa and Pc partitions of the rows and the columns of A and B [§, the 
eigenvalue equation can then be rewritten as 



det[u — S — ^{^)) = minor{ui — S) minor (E,{u!)) 

Pr,Pc 



^ det{K^^y) 



(5.15) 



In the second line of Eq. (|5.15p we considered the two terms in the minor ex- 
pansion that have the maximum and the minimum degree in a;, respectively 
Ns and —Ng. Thus, assuming a completely general structure for the if^^^ 
terms, Eq. ()5.15p is a polynomial equation of degree Ng + NdNg. Conse- 
quently the introduction of a frequency-dependent kernel yields, in general, 
more solutions then the single electron transitions (Ng), although larger then 



■^We recall that a minor of a matrix A is the determinant of a submatrix M obtained 
from erasing a fixed number n of columns and rows. The terms n = 0, i.e. the determinant 
of the matrix M = A, is considered too. 
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the correct number of poles, Ns + Nd- However, in our case, the particular 
structure of the matrices K^^\j^ = Cij^^C^j. ^ ensures that the determinant of 
any but the one-dimensional sub-block of K^'^'> is zero. This means that the 
second term on the r.h.s. of Eq. ()5.15|) is zero and in the minor expansion 
only Nd terms of degree —1 survives, from which it follows that the total 
degree of the polynomial expression det{uj — S — 2(a;)) is Ng + Nd- 

In the notation of Eq. (|4.12p ei = A\^eh- By plugging the eigenvectors 
and eigenvalues of Eq. (|5.8|) in Eq. ()4.12|) we see immediately that the Nd 
DEs will appear as poles of ^(oj). 

Feynman diagrams reteined in the second random phase 
approximation 

In order to create a common language between the sRPA and the DBSE ap- 
proaches we start by noticing that, within TDA, the kernel of sRPA contains 
all Feynman diagrams up to second order. The 16 second order diagrams 
included are represented in Fig. (|5.ip . 

To understand why DEs are described within this approximation we focus 
on diagram (a) of Fig. ()5.ip . drawn for a specific time ordering. The diagram 
describes a physical process where the electron-hole pair created at time t 
emits a photon that generates another electron- hole pair at time ti. The 
second e-h pair is annihilated at time t2- Therefore this Feynman diagram 
is describing the coupling between a single- and a double-excitation. 

Second random phase approximation, correlation, and 
TDA: a closed end 

In extended systems the dressing up of bare particles induced by correla- 
tion effects is mediated by collective charge oscillations, i.e. by plasmons. 
Therefore a coherent approach to DEs in correlated materials should also 
describe the interaction with plasmons. The key problem in the description 
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(b) e-h exchange diagrams 




(c) one particle exchange diagrams 



Figure 5.1: Basic Feynman diagrams included in the second-RPA approach [531 
154) beyond the standard TDHF. The time flows from left to right respecting the 
Tamm-DancofF approximation. The sRPA approach, when the TDA is relaxed, 
includes other 16 basic diagrams obtained by inverting the direction of all GF. The 
complete set of diagrams is obtained by iterating the Dyson equation. 



of plasmons is the possible breakdown of the TDA, as it oceurs, for exam- 
ple, in nano-structures |56| . Indeed, within the TDA neutral excitations are 
described as packets of electron-hole pairs propagating only forward in time, 
and, therefore, charge oscillations (plasmons) cannot be captured. 

sRPA can, in principle, describe plasmons by going beyond the TDA. 
However, as a matter of fact, the complexity of the method imposes to 
retain only a few terms beyond TDA. Indeed, the sRPA, given by Eq. (j5.8|) . is 
equivalent to a Dyson equation for the response function that can be analysed 
by using the diagrammatic technique. It results that while kernel diagrams 
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(a) kernel diagrams 




(b) self-energy diagrams 



Figure 5.2: Second order Feynman diagrams relevant to the description of col- 
lective excitations. Time flows from left to right. While kernel diagrams (a) are 
included in the sRPA, trough the iteration of the Dyson equation, self-energy ones 
(b) are not. This inconsistency prevent the sRPA to work in a correlated regime. 



(see Fig. (|5.2p . panel (a)) arc included in tlie sRPA, self-energy diagrams (see 
Fig. (|5.2p . panel (b)) are not. It has been shown that, starting from the HF 
approach, including only the kernel diagrams yields an incorrect description 
of the excitation energies [52 • In a recent paper by Gambacurta et al. [15], 
studying the spectrum of Sodium clusters, this problem is discussed and 
identified as lack of ground state correlation. The same problem is identified 
by Huix-RotUant and Casida [32] . 

This is one of the major reasons why the sRPA approach is not very 
popular in the condensed matter field. Approaches like the Algebraic Dia- 
grammatic Construction (ADC) are preferred |361[55]. However in the ADC 
approach kernel and self-energy diagrams beyond TDA are included only 
up to finite orderfl While this is a reasonable approach for small systems, 
it is expected to fail in extended correlated ones. In extended system any 

•^Another limit of the ADC scheme, from our point of view, is that this explicitely 
includes all diagrams related by particle exchange and Pauli exclusion principle. For 
double excitations this forces to perform a matrix diagonalization in the space of 2p — 
2/i, when the electron hole interaction among virtal particles is considered (ADC(2)— x 
approximation 1361 ) 
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order diagram in the bare interaction is relevant and kernel and self-energy 
diagrams must be included up to an infinite order. 

5.2 The dynamical Bethe— Salpeter equation 
(step II) 

It is now clear that a well defined approach to the description of DEs must 
be NC, i.e., it must not introduce spurious non-physical solutions. At the 
same time it must include diagrams up to infinite order and beyond the 
TDA in order to describe screening effect and collective excitations. The 
BSE approach is an alternative scheme which includes the infinite series of 
both kernel and self-energy diagrams, thus providing a suitable approach to 
achieve both goals in a coherent manner. 

However in the DBSE presented in Ch. 2] the kernel S2(w) includes 
the frequency dependency of kernel diagrams only, whereas the static self- 
energy effects are included as a rigid shift of the QP eigen-energics. It 
has been already shown that, at linear order, dynamical effects have to be 
included in both the kernel and the self-energy [Sn|. Following this in- 
put, in order to construct a consistent dynamical approximation, we im- 
prove the DBSE including dynamical self-energy effects in a term 
to be added to the kernel ^^(w). This extra term originates from the 
Lo{u!) — —i J duj' / {2tt)G{ll!' + uj /2)G{u!' — u! /2) as described in the following. 
We use the Dyson equation for the GF written in the form: 

= g-^ - - J^diuj) = G;^ - (5.16) 
where we separated the Self-Energy in its static Ss and dynamic parts. 
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and we defined Gg — g — Sg- In this way we can write: 
Lo(w) = Ls{uj) 

/J I 
-^Gsiuj' + uj/2)Gs{uj' - uj/2)T.d{uj' - uj/2)G{uj' - uj/2) 
2-K 

/I ! 
-^Gsiuj' + a;/2)I]d(tj' + cj/2)G{uj' + uj/2)Gs{uj' - uj/2) 
2n 

/I f 
-^Gsiuj' + uj/2)T.d{uj' + uj/2)G{J + uj/2) 
2-K 

G,{J - cj/2)Sd(w' - w/2)G(w' - w/2), (5.17) 

where = —iGsGg- Using the same trick adopted for the kernel, we mul- 
tiply the second, third, and fomth term on the right-hand side of Eq. (|5.17|) 
by LsL~^ from the left and by Lq^Lq from the right, and wc obtain 

Lo(w) = i.(w) + L,{iv)Ei{io)Lo{oj) (5.18) 

with 

[Gsiuj' + uj/2)Gsiuj' - w/2)Erf(w' - uj/2)G{uj' - lj/2) 
+ Gsiu}' + w/2)Ed(w' +uj/2)G{uj' + w/2)Gg(w' - w/2) 
+ G,{uj' + tj/2)Sd(w' + w/2)G(w' + Lo/2) 

Gs{uj' - w/2)I]d(w' - uj/2)G{uj' - w/2)] Lq^{uj). 

For the description of DEs we set to zero the last term on the right-hand 
side of Eq. (|5.17p as it describes a process where six Green-function lines 
appear in the same moment, so a triple excitation. We thus obtain a total 
kernel S''(a;) — Sf (oj) -f- "BJ^iio) to be inserted in Eq. (|4.17p . which wc write 
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here again for clarity 



(5.19) 



The zero order term is now called Lg[uj) to underline that it includes only 
the static effects of the self energy. 



5^(w) = — -- / dLu'duj"dLo"'Laiuj,io',Lu") 

{2TTy J 

E{lj, lj\lu"')L{uj, w'", lo")L~\uj). (5.21) 



The complexity of the original Eq. (|4.13p is thus transferred in the structure 
of the DBSE kernel E!d(w). We can, however, simplify the dependence on 
G and L in '^di^) by starting from its linear limit where G{uj) ~ Gs{uj) in 
Eq. ([OO)) and L{uj,uj',u") ~ Lq{uj,uj' ,uj") ~ Ls{uJ,uj' ,uj") in Eq. (jOTj) . 
This limit is fully justified in the DBSE by the fact that it accounts for 
the simultaneous evolution of two c-h pairs, which represent the dominant 
channel in the description of DEs. Accordingly it is crucial that the static 
part of the self-energy is treated in a separate way. In contrast to Tfi in fact, 
the static part of the self energy cannot be treated through a linearized 
kernel, as this would lead to numerical instabilities |60| . 

We need now to approximate the unknown quantities E^, Yid{<^) and 
S(a;, a;', oj"). One can verify that the DBSE is equivalent to sRPA if one 
chooses Tig = T,hf, "^d = '^2{^), thus including all second order Feynman 
diagrams of sRPA, and S(ti;, w', uj") ~ dk'^B^^i^^)) ^ gi^^j-^jj^g from this obser- 
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vation in the next section we derive a diagrammatic number conserving rule 
which we then use to construct an approximation able to properly describe 
correlated system. To do this we will consider the GW approximation to 
the Self-Energy evaluated at the QwP eigen-energies and we will derive the 
dynamical part and the kernel E!(a;, w', w") starting from the screened 
Coulomb interaction, in order to include the static GW-BSE scheme in the 
— > limit. The derivation will be carried on within TDA in order to keep 
the discussion as simple as possible. 

5.3 A number— conserving kernel for correlated 
systems 

The DBSE equation provides a powerful starting point to tackle the double- 
excitation problem, as the diagrammatic approach makes possible to intro- 
duce different levels of approximation that overcome the limits of the sRPA. 
We achieve this by following two essential steps: i) we use the sRPA to 
create a close link between the diagrams introduced in the DBSE kernel, 
the particle indistinguishability and the Pauli exclusion principle; ii) we use 
this link to define a number-conserving correlated kernel starting from the 
standard GW approximation. 

The diagrammatic number conserving rule 

By taking into account the 16 diagrams of Fig. (|5.ip the DBSE (and conse- 
quently the sRPA) correctly describes particle indistinguishability and Pauli 
exclusion principle. Here we illustrate how this can be deduced from the 
inspection of the Feynman diagrams. The 16 diagrams describe processes 
in which a DE appears from a photon emitted either from the electron or 
from the hole and then absorbed back (these two possibilities are the first 
two terms in the definition of C(,;j) („g)(mp), see Eq. (|5.1ip ). Therefore each 
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(a) e— h pairs exchange 



(b) Single particles exchange 



(c) Half diagrams recombination 



Figure 5.3: The DNCR in practice. We take as reference the two time orderings 
of the kernel bubble diagram, corresponding to the two first diagrams of Fig. (|5.ip . 
The general procedure to get a NC kernel is to split each initial diagram in two 
half-diagrams. Then these half diagrams must be connected by exchanging in all 
possible ways all e-h pairs and all single particles. This produces a new group 
of diagrams that must be processed using the same procedure. When no new 
diagrams appear the resuhing kernel is NC. 



double-excitation process can start and end in two ways so that there are 
4 possible processes, which are the four bubble diagrams of Fig. ()5.ip . The 
other 12 diagrams reflect the particles indistinguishability that imposes the 
electron lines, as well as the hole lines, to be interchangeable among them- 
selves. 

Following Fig. (j5.3p we can derive a graphical rule that any approxima- 
tion has to respect in order to be NC , this is the proposed Diagrammatic 
Number Conserving Rule [DNCR). First we consider an initial group of 
diagrams, chosen in such a way to describe the relevant physics we want to 
introduce in the theory (like plasmons and excitons) . Then we split each dia- 
gram in two parts that, connected in all possible ways obtained by imposing 
particles exchange, lead to a new group of diagrams. When the same pro- 
cedure applied to the resulting diagrams does not lead to any new diagram. 
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then the approximation is, by definition, NC. As an illustration, the DNCR 
can be applied to the sRPA diagrams, shown in Fig. (|5.ip . It can be shown 
that all sRPA diagrams can be obtained from the first two by applying the 
proposed DNCR. 
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Figure 5.4: Basic RPA diagram used as a starting point for the correlated kernel. 
All other diagrams are obtained by applying the DNCR, as discussed in the text. 

The diagrammatic number conserving rule applied to the 
Bethe-Salpeter equation 

A crucial consequence of the DNCR is that, as exemplified in Fig. (|5.3p . a 
NC kernel must include all kind of diagrams. Therefore, whatever initial 
approximation is chosen the repeated application of the DNCR will create a 
balanced mixture of diagrams in order to respect particle indistinguishability. 
If the DNCR is not respected by selecting only a class of diagrams, then 
spurious solutions are expected to appear. This is the case of the kernel 
proposed in Ch. IHthat was obtained from the standard S ~ ly (ui) by simply 
relaxing the static approximation for W. This kernel introduces an infinite 
series of RPA diagrams only in the interaction W, neglecting all consequent 
diagrams imposed by the DNCR. As a consequence spurious poles in the 
polarizability are found as predicted by the DNCR. 

Nevertheless the kernel proposed in the previous chapter describes the 
interaction with plasmons, which is a desirable property which we want to 
retain, at the same time forcing the kernel to be NC. However, before ap- 
plying the DNCR, we have to note that the W{uj) propagator describes the 
evolution of charge oscillations, composed by renormalized packets of e-h 
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Figure 5.5: Build up of correlated Feynman diagrams connecting two Feynman 
diagrams. The effect of exchange among two RPA excitations is shown using dotted 
lines. 



pairs. This clearly makes a distinction between the e-h pairs embodied in 
W (u) and the real e-h pairs created by the scattering process leading to 
the breakdown of the particle indistinguishability. A better starting point 
is instead the basic diagram showed in Fig. 15.41 where all e-h pairs are cor- 
rectly renormalized. In this diagram the filled bubble and the filled rectangle 
represent the RPA response function t^^^^ (o;). By introducing the Lehman 
representation for x^^^ in the same notation of the previous chapter we can 
write 

xl%^i^)=T. (5-22) 

We will call the poles of x^^^ RPA excitations. Note that these capture 
the physics of the plasmonic oscillations. 

The DNCR imposes to consider all possible diagrams obtained from 
Fig. 15.41 by exchanging the basic excitation propagators. The key point 
here is to rotate from the independent e-h pairs to the RPA basis, where 
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e-h pairs are replaced by the RPA excitations. Therefore we proceed by 
spUtting the RPA propagators, using Eq. (|5.22p . as sketched in Fig. (|5.5p . 
Then we consider ah diagrams where the RPA excitations are exchanged. 

Mathematically the procedure sketched in Fig. ()5.5|) corresponds to rotate 
in the RPA excitation space the residuals and poles of Eq. ()5.10p . Each 
term in the rotated counterpart of Eq. (jS.lip will correspond to a possible 
connection induced by the DNCR: 

nq.mp ^ 

(5.23) 

As a consequence the correlated version of Eq. ()5.10p will look like 

(^B.PA \f<RPA 

The symbol {vi o V2] in Eq. (|5.23p imposes the invariance of the correlated 
kernel under exchange of RPA excitations. Consequently the kernel S^p^ 
is by definition invariant under exchange of two RPA excitations. However 
RPA excitations are bosons so that Pauli exclusion principle is not taken 
into account and the obtained kernel is not fully NC. To fix this problem it 
is sufficient to impose the condition vi 7^ in Eq. (|5.24p . 

The DBSE obtained by using the ^^j^pj^ kernel includes all self— energy 
terms obtained from = GW and S(w, w', u") ~ (5(cj — cj')VK(cj' — w"). In 
addition extra terms appear in order to fulfil the NC condition. Interest- 
ingly Spp^, within TDA, also embodies the full frequency dependent term 
G SW 1 5G(uj) which is usually neglected in the standard BSE approach. In 
the present case these second order diagrams in W are indeed needed to 
correctly account for the particle indistinguishability. The resulting kernel, 
whose diagrammatic expression is sketched in Fig. (|5.6p . does have the right 
mathematical structure by construction, so that no spurious solutions are 
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Figure 5.6: The final correlated DBSE kernel. The filled regions represent the 
propagation of RPA excitations. 



present. 
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5.4 Numerical results on model molecular sys- 
tems 

In the following wc will illustrate various conceptual and technical aspects 
of our approach using two benchmark model systems, based on the Cg,H2 
and the Ca^H^ molecules. These unsaturated hydrocarbon chains are often 
chosen as benchmark systems to test theoretical methods aimed to describe 
double excitations. By calculating the polarizability of these systems we 
will show: i) the role played by subgroups of diagrams in the description 
of double excitations; ii) the fact that the number conserving rule not only 
applies to the total number of poles, but also to the number of optically 
active poles; iii) the absence of spurious double excitation peaks that appear 
that appear in approaches [33] that violate the NC rule. 

The calculations have been performed using the YAMBO code [ST], where 
we implemented sRPA for closed-shell systems, within the TDA. Further- 
more we approximate both QP and HF wave-functions with KS-LDA wave- 
functions. 

sRPA produces results similar to the GW-BSE approach or to the DBSE 
when only "bubble diagrams" (first row of Fig. 15. ip or bubble diagrams and 
"eh exchange diagrams" (second row of Fig. 15. respectively, are selected. 
Therefore this implementation allows us to explore the performances of the 
various approaches by selecting specific subgroups of diagrams. 

We first performed a ground-state calculation with the ABINIT code [52] , 
within DFT/LDA, with an energy plane-wave cut-off of 20 Hartree and a 
super-cell of 25 x 25 x 40 Bohr for the C%H2 (a linear molecule w 21 Bohr 
long) and a smaller super-cell of 25 x 25 x 15 Bohr for the C^Hq (the 
molecule extends for sa 10 Bohr both in the x and y directions). Then we 
performed excited-state calculations in the basis-set of KS-states, consider- 
ing only the states from HOMO-3 to LUMO+3. In this way our systems can 
be mapped into an eight level model with 16 single and 240 double excita- 
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5 6 7 8 9 10 11 12 13 14 5 6 7 8 9 10 11 12 13 14 



(a) c,H^ energy [eV] (b) c,H^ 



Figure 5.7: Second Random Phase Approximation spectra. For botli model sys- 
tems tiie frequency-dependent Icernel produces extra peaks (red line) which cannot 
be described by a static kernel. The black thin dashed line is the Independent- 
Particle spectrum. The inset is present here as reference to detect spurious peaks 
in the insets of Fig. 15.81 and 15.91 Alia spectra in this and in the following figures 
include an artificial broadening due to the use o an imaginary factor ir] = O.OSeV 
in the Green's function denominator. 



tions. All the C%}i2 eigenvalues are doubly degenerate due to the symmetry 
of the molecule. 

In the description of double excitations the kernel frequency dependence 
becomes crucial when one or more poles fall in the absorption spectrum 
energy range. In this case the static approximation fails, and extra peaks 
appear. In order to artificially simulate this situation in our systems wc use 
HF eigenvalues to construct ip: while the kernel is built with KS-LDA ones. 
This choice gives us the possibility to investigate more physical situations 
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which could arise for correlated materials. 

The results of these calculations are plotted in Fig. 15.71 For both systems 
at the HF Independent-Particle (IP) level there is a clear peak, which falls 
close to 9 eV for the CsH2, and close to 12 eV for C^Hq. As expected 
the kernel constructed with KS-LDA eigenvalues has poles in these energy 
ranges, so that extra peaks appear in the spectrum. The effect is visible in 
both model systems: in Cg,H2 the main peak is essentially split in two (see 
Fig. 15.71 panel (a)); for C^Hq several extra peaks appear as shown in the 
inset of Fig. 15.71 panel (6). 

Wc will now explore the role played by the various subgroups of diagrams, 
namely (a) the "bubble diagrams" (first row of Fig. 15. 8^ : (b) the "e/i exchange 
diagrams" (second row of Fig. 15. 8|) . which are obtained from the bubble 
diagrams via eh exchange; (c) the "particle exchange diagrams" (third and 
fourth row of Fig. 15. 8p . which are obtained from the bubble diagrams via 
single-particle exchange. 

Fig. 15.81 shows the spectra obtained taking into account, beyond the 
TDHF scheme, only selected families of diagrams. By selecting only dia- 
grams of type (b) or (c) the spectra are not positive defined. This unphys- 
ical property can be understood by noticing that the frequency-dependent 
kernel constructed from diagrams (b) and (c) does not have the mathemat- 
ical structure of Eq. (|5.10|) . On the contrary the kernel constructed from 
the bubble diagrams (a) is positive, though particle indistinguishability and 
Pauli exclusion principle arc not respected as illustrated in previous sections. 

Indeed the spectra constructed from bubble diagrams is positive defined, 
though spurious peaks appear: in Cg,H2 (Fig. 15.8] panel (a)) one has three 
peaks at around 12 eV, and in the C^^H^ many peaks appear (see Fig. 15.81 
panel (6), the left inset) which arc not present in the full sRPA spectra. 

Fig. 15.91 shows the spectra constructed taking into account both subsets 
of diagrams (a) and (b) or (a) and (c) together. The spectra are positive 
defined. However, only the former combination gives the right number of 
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(a) CgH^ (b) C^H^ 



Figure 5.8: The spectra obtained selecting only specific subgroup of diagrams. 
By selecting only "eh exchange" (green dashed line) or "one particle exchange" 
(blue dots) an unphysical (negative) polarizability is observed. Only the spectra 
obtained with the kernel constructed using the "bubble diagrams" (red line) is 
positive defined. However when only "bubble diagrams" are used, as proposed by 
Romaniello et. al [33], spurious peaks appear. These peaks do not appear in the 
spectra obtained from the full sRPA kernel (see Fig. 15. 7p . 

peaks (i.e. the same number of the full sRPA spectrum) whereas the latter 
produces spurious poles. 

In this perspective it is interesting to compare the two cases. In the 
6*8^^2 model the subset of diagrams (a) and (b) (green dashed line) gives 
a spectrum which is very close to the full sRPA spectrum of Fig. 15.71 (red 
line) both in the structure and number of poles. Diagrams of kind (c) are, 
instead, negligible. 

In the C^Hq model, on the contrary, diagrams of kind (c) play an im- 
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Figure 5.9: Spectra obtained witii tiie kernel constructed using tiie "bubble+eh 
excliange" diagrams (green dashed line) and the "bubble+eh exchange" diagrams 
(blue line). Both spectra are positive defined but only the combination "bubbles 
+ eh exchange" gives the same number of poles of the full sRPA spectra (see 
Fig. I5.7|l . On the contrary the combination "bubbles + one particle exchange" 
gives spurious solutions. 



portant role: they shift the peak of the "bubbles" polarizability towards the 
results obtained with the sRPA kernel. Diagrams of kind (b), instead, in 
this case have a negligible effect on the position of the peaks. However the 
choice (a) + (c) gives several spurious poles (see Fig. I5.9l panel (6), blue line 
in the right inset) and, as for CsH2, only the combination (a) and (b) yields 
the correct number of poles (see Fig. 15.91 panel (6), green dashed line in the 
left). 

The sum of diagrams (a) and (b) describes eh pairs as indistinguishable 
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bosons, whereas the sum of diagrams (a) and (c) does not correspond to any 
defined statistic. We can then conclude tliat diagrams of kind (c) are mean- 
ingful only if added to the other two classes of diagrams in order to describe 
particle indistinguishability. However, the spectrum obtained combining the 
diagrams (a) and (b) has indeed the same number of peaks of the spectrum 
obtained using the complete kernel, thus supporting our recipe to construct 
a correlated kernel discarding the subset of diagrams (cjf|. Another conclu- 
sion we draw from these results is that our approach, by respecting the NC 
rule, ensures that the theory produces not only the correct total number of 
poles, but also the correct total number of optically active (and optically not 
active) poles. 



■^The subset of diagrams (c) could be included in a correlated kernel only at the price 
of a direct diagonalization in the space 2p — 2h as in the ADC(2)-x scheme 
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5.5 Conclusions 

In this part of the thesis we presented a method to include double excitations 
in a consistent manner within the GW+BSE approach. The main idea has 
been to correct the standard BSE kernel in order to go beyond the static 
approximation fulfilling the Number Conserving (NC) condition. The result- 
ing scheme keeps all the advantages of the Many Body approach, that is the 
ability to describe extended and correlated materials in a consistent man- 
ner, without producing spurious excitations, provided the cndition vi ^ vi 
in Eq. (|5.24p is imposed. This is not an exact condition of the kernel but is 
sufficient to give the correct numbering of poles. 

The NC condition results from an inspection of the similarities and the 
differences between the BSE scheme, designed for solids, and the sRPA ap- 
proach, designed for isolated systems. The main character in the first is the 
screening, while in the second it is the role of exchange. As pointed out in 
the very recent work by Huix-RotUant and Casida [32], there is a great in- 
terest in this direction in order to develop approximations at the nanoscale 
interface between molecules and solids. As opposed to other works, how- 
ever, we do not directly consider all exchange diagrams related to the RPA 
screening re-summation, because we believe that such an approach would be 
impractical, especially for nanostructured materials. Instead, our method is 
aimed to capture the main feature related to the exchange principle without 
requiring matrix diagonalizztion in the space of double excitations. 



Part III 



Carbon nanotubes in 
magnetic fields 



Chapter 6 



The Aharonov Bohm 
effect in carbon nanotubes 

The Aharonov-Bohm (AB) effect [S3] is a purely quantum mechanical effect 
which does not have a counterpart in classical mechanics. A magnetic field B 
confined in a closed region of space alter the kinematics of charged classical 
particles only if they move inside this region. Electron dynamics, instead, 
governed by the Schrodinger equation, is influenced even if the particles 
move on paths that enclose the region where the magnetic field is confined, 
where the Lorentz force is strictly zero. If this closed region is the inner 
part of a nano-tube, electrons traveling around the cylinder are expected 
to manifest a shift of their phase. The mathematical interpretation of this 
effect is connected with the definition of the vector potential, which, in the 
case of confined magnetic fields, cannot be nullified everywhere. 

This extraordinary effect, first predicted by Aharonov and Bohm [S3] (AB) 
in 1960, was interpreted as a proof of the reality of the electromagnetic poten- 
tials. The idea that electrons could be affected by electromagnetic potentials 
without being in contact with the fields was skeptically received by the sci- 
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entific community. At the same time the AB paper spawned a flourishing of 
experiments and extension of the original idea. The first experiment aimed 
at proving (or disproving) the AB effect revealed a perfect agreement with 
the theoretical predictions [64]. Nevertheless only some years later, in 1986, 
the experiment which can be considered as a definitive proof of the correct 
interpretation of the AB effect was realized. Tonomura et al. [SS], using 
superconducting niobium cladding, were in fact able to completely exclude 
the possibility of stray fields as alternative explanation of the predicted and 
observed AB oscillations. 

Nowadays the AB effect can be used in a wide range of experiments, from 
the investigation of the properties of mesoscopic normal conductors to the 
experiments designed to reveal the structure of flux lines in superconductors. 
Growing interest is emerging in the field of nanostructured materials. One 
of the most well-known case is given by Carbon NanoTubes (CNTs) that, if 
immersed in a uniform magnetic field aligned with the tube axis, have been 
predicted to show peculiar oscillations of the electronic gap. These oscilla- 
tions are characterized by a period given by the magnetic fiux quantum h/ e 
and are commonly interpreted as caused by the change in the wave functions 
of the electrons localized on the tube surface induced by the Aharonov-Bohm 
effect. 

The first experiment carried on CNTs, in 1999, described the oscillations 
in the electronic conductivity |66| . but with period of h/2e. This devia- 
tion from the predicted AB oscillation period has been explained in terms 
of the weak localization effect [ST] induced by defects and dislocation by 
Al'tshuter, Aronov, Spivak [68] (AAS effect). Only in 2004 a clear proof 
of the existence of AB oscillations with an h/e period have been given by 
Coskun et coll. |69j by measuring the conductance oscillations in quantum 
dots. The dots were built using concentric Multi-Wall (MW) CNTs of dif- 
ferent radii, short enough to prevent the appearance of weak localization. 
In the same year Zaric et al. [70] observed modulation in the optical gap of 
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(a) AB picture (b) AB experiment 



Figure 6.1: Representation of the AB effect. In the left panel, panel (a), in 
green the vector potential generated from a solenoid in the region B = in the 
symmetric gauge. From Ref. [71]. On the the right the ideal experimental setup 
of the AB effect. From Ref. [03]. 

pure Single Wall (SW) CNTs with oscillations of h/e period. 



6.1 What is the Aharanov— Bohm effect? 

The AB effect was introduced in |63| by considering the interference exper- 
iment described in Fig. (|6.1I &). An (ideally infinite) solenoid generates a 
magnetic field only inside the solenoid itself, hi contrast to the magnetic 
field, the vector potential A, which satisfy the condition V x A = Bq inside 
the solenoid, will not be zero outside. Indeed in the symmetric gauge the 
total vector potential written in cylindrical coordinates, A — (0,0, A^) is 

= l-Bor for r < tq, 

I ,2 (6-1) 

A^^-Boy forr>ro. 
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where ro is the radius of the solenoid. In Fig. (|6.1l a). the vector poten- 
tial around the solenoid is represented by the green circumferences. In real 
experiments, as the solenoid is not infinite, the return magnetic field must 
properly set-up in order to avoid the regions where the electrons are per- 
mitted. Electrons are injected in the experiment from a point outside and 
far from the solenoid (see Fig. 16.11 6) 

If electrons were classical particles their equation of motion would be 
completely determined by the magnetic field only and in the experiment 
of Fig. (|6.1l a) they would not be affected at all by the presence of the 
solenoid. Electrons, instead, are quantum particles and arc governed by 
the Schrodingcr equation, where the potentials do enter. In |63j Aharonov 
and Bohm demonstrated that the electrons feel the presence of the solenoid 
by acquiring a phase shift between the two paths of Fig. (|6.1l a). This phase 
can be measured as change in the interference pattern on the screen. The 
phase shift, >S'(x), between the two paths can be expressed in terms of the 
vector potential as 



and, computing the integral along the closed line obtained from the path of 
the two wave-packets, we obtain S'(x) = e//i $ with $ the magnetic flux. 
The electronic wave-function which describes the image on the screen is then 
*' = *o e-'*/*«, with *o the wave-function when the experiment is carried 
out without the solenoid and ~ h/e the flux quantum. 

A static version of the AB effect also exists. If we consider a free electron 
on a ring which encloses a solenoid we can compute its eigenvalues as |72) 




(6.2) 




(6.3) 



where is the canonical angular momentum of the electron, r is the radius of 
the ring and $ the magnetic flux through the ring. We see that the presence 
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of a magnetic flux modifies tlie eigenvalues splitting their ±lz degeneracy. 

The quantum mechanical nature of the AB effect is clearly shown by its 
proportionality to the magnetic flux $0i which goes to zero in the classical 
limit. 

Interpretation of the Aharonov— Bohm effect 

When they first proposed the existence of the effect, Aharonov and Bohm 
claimed that it was a proof of the reality of the electromagnetic potentials. 
Their paper generated an intense debate in the scientific community which 
is not yet terminated. Here we offer some considerations which are inspired 
by the review by Peshkin and Tonomura j72j on the AB and in particular 
on the "central role of the quantized angular momentum". 

Consider the static AB effect previously described. We will show in 
the following that, if the AB effect did not exist, than we would obtain, 
as a result, that the eigenstate of the system depend on the history of the 
system. This is in sharp contrast with the foundations of quantum mechanics 
which states that the Hamiltonian of a system is, at any time, a well defined 
operator with unique eigenvalues and eigenvectors. 

The hypothesis that the AB does not exist means that and we can com- 
pute the electronic eigenstates in a (multi-connected) region of space from 
the sole knowledge of the magnetic field in that region. This means, in our 
static example, on the ring and nearby it where the electronic wave-function 
is different from zero. 

Suppose that at an initial time t = to there is no current flowing through 
the solenoid and one electron is in a steady state of the ring, Ei = 2^72 
with total kinetic angular momentum in the z-direction = Mz- Then 
we turn on the current and, during the transient, a time dependent electro- 
magnetic field is generated by the solenoid. This field will generate a torque 
on the electron AK^ = fi{e^/h) (details of the calculation can be found 
in [72j ) so that the total angular momentum of the electron, which must 
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then be conserved, becomes = ~ e^/h). 

Let's now consider a second possibility, where the electron on the ring 
originates from an eh pair created long after the electromagnetic wave has 
been dissipated. In this case the electromagnetic field is zero on the ring and 
nearby. Then the steady states on the ring, accordingly to our hypothesis 
that the AB effect does not exist, have eigenvalues Ei = ^^^^ indepen- 
dently of the presence of the solenoid. 

The natural conclusion would, then, be that the electronic eigenstates 
depend on the history which is clearly in contradiction with the principles of 
quantum mechanics. As a consequence the AB effect must exist in order to 
ensure that the eigenstate in presence of a magnetic field, possesses angular 
momentum = fi{lz ~ e^/h). Thus we can interpret the AB effect as 
a witness of the previous switching on of the solenoid, which modified the 
space around itself. 

Persistent Currents 

The existence of the AB effect is strongly related to the quantization of the 
angular momentum in quantum mechanics [72j : to be precise, in the Hamil- 
tonian formalism the canonical angular momentum L = r x p is quantized. 
When a vector potential A exists, the mechanical angular momentum of 
the electron is K = r x (p — A) which is, in general, non quantized. This 
observation has an important consequence: the existence of Persistent Cur- 
rents (PCs) in quantum mechanics generated when clccrons move in some 
particular topologies, like rings or cylindrical shaped objects. 

In classical mechanics if we move a metallic ring in a region of space where 
a magnetic field is present the change in the magnetic fiux induces a transient 
voltage and a current, that will eventually disappear due to the existence of 
dissipation mechanisms. The appearance of an angular current proportional 
to the vector potential j^^^ = -\-e^/rn? Aj^Pp is in fact counterbalanced by 
an opposite current j'^P'' = ~e/ra Re\9* p 5*] that relaxes the system to the 
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lower energy configuration witli j = 0. Tlic current flowing in tlie angular 
direction is proportional to the angular momentum divided by the radius of 
the ring. However the term r x p can assume only integer values and the 
total current can be nullified only when r x A is an integer, while for any 
other value PCs exist. If we consider a constant magnetic field then, in the 
symmetric gauge, 

r = (6-4) 

can be nullified only when $/$o is an integer. Here p is the electronic density 
and Rring the radius of the ring where electrons are trapped. The existence of 
a periodic lattice partially breaks the quantization of the anguar momentum, 
however in mesoscopic rings (or cylinders) the quantization is almost exact 
and even at the nanoscale the argument is correct in first approximation. 
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The possibility of observing PCs in non superconducting mesoscopic 
metallic rings was first proposed by Biittiker et al. in f983 [73], while ex- 
perimental confirmation was reported in f990 by a research group at Bell 
Laboratories [75]. Fig. (|6.2|) provides a schematic representation of PCs in a 
mesoscopic ring. PCs have been predicted to exist in nanostructured mate- 
rials, like carbon nanotubes |76| . However PCs have been measured only in 
2009 in gold and aluminum rings [77] [75]. PCs in fact are easily destroyed 
increasing the temperature of the system by the smearing on the electronic 
occupations. In particular if the smearing is greater than the difference in 
energy between states with opposite angular momentum the currents vanish; 
for this reason PCs do not exist in macroscopic objects. 

PCs have never been measured in CNTs. We will discuss this subject 
more in details in the next chapter. 

The Al'tshuler, Aronov and Spikav effect 

In the next section we will describe how the AB effect influences the elec- 
tronic properties of CNTs. We will provide some experimental evidences, 
like the resistivity oscillations observed in MW-CNTs with a period which 
is half the AB period. To understand this result we need to introduce here 
the concept of Weak Localization (WL), which was first proposed in the 
1979 [67]. WL is often seen as a precursor of strong localization in dis- 
ordered materials and the detailed derivation of the phenomena is rather 
intricate while its phenomcnological interpretation is quite intuitive. 

In the classical theory of transport phenomena the total probability for 
a particle to transfer from point P to point Q (Fig. (|6.3|) ) is the sum of 
probabilities of such a transfer over all possible trajectories. In quantum 
mechanics this result corresponds to neglect the interference of scattered 
electrons propagating along different paths and having approximately ran- 
dom phases under the quasi-classical condition X << I, with / the length of 
the propagation path and A the De Broglie wave-lenght of the wave-packet. 
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Figure 6.3: Different types of quasi-classical particle trajectories connecting P 
and Q. Point O is a self-crossing point of two possible trajectory. From Ref. [75] 

There is, however, a specific class of trajectories, namely, self-crossing tra- 
jectories (trajectory 2 in Fig. (|6.3p ) for which the wave interference turns 
out to be essential. Indeed two waves propagating along such trajectories in 
two opposite directions (conjugated waves) accumulate the same phase dif- 
ference. Therefore the contribution of these trajectories to the probability 
of coming to the same point (point O in Fig. (j6.3p ) will be 

|Ai + = + |^2|2 + 2Re[AlA2] (6.5) 

which is twice the sum of the squared amplitude moduli. A higher probability 
of returning back to point O means a lower probability of transfer from point 
P to point Q. Thus (weak) localization is favourite and, hence, results in an 
increase of the resistivity. 

If the sample is placed in a magnetic field then the probability amplitudes 
of completing the loop on contour 2 of Fig. (|6.3p acquire an additional phase 

Ai -> Aie^'""^ A2 Aae'^"^ (6.6) 
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and then the phase difference wili be ^ip/{2'n) = 2$/$o- Indeed the com- 
bined WL and AB effects predict the existence of resistivity oscillations with 
a period of $o/2 and is known as AAS effect from the names of Al'tshulcr, 
Aronov and Spikav who first proposed its existence in 1980 [68) . 



6.2 An introduction to carbon nanotubes 

A CNT is a honeycomb lattice rolled into a hollow cylinder with nano-metric 
diameter and /zm length. CNTs were discovered and first characterized in 
1991 by lijima from NEC laboratories (Japan) [1]. The first CNTs discovered 
were made of several concentric cylindrical-like shells regularly spaced by an 
amount of about 3.4 A as in conventional graphite materials. These Multi- 
Wall CNTs (MWCNTs) were first synthesized with diameters ranging from a 
few nanometers to several hundred nanometers for the inner and outer shells, 
respectively. As for the length, MWCNTs extending over several microns are 
currently synthesized. Shortly after the discovery of MWCNTs, Single- Wall 
CNTs (SWCNTs) were synthesized in abundance using arc-discharge meth- 
ods with transition-metal catalysts [HH [H2] • These tubes have quite small 
and uniform diameter, on the order of Inm = lO^^m. This unprecedent- 
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edly small diameter, combined with the crystalline perfection of the atomic 
network, explains why these objects were quickly considered as the ultimate 
carbon-based ID systems. Crystalline ropes or bundles of SWNTs, with 
each rope containing tens to hundreds of tubes of similar diameter, closely 
packed in a hexagonal configuration, have also been synthesized using a laser 
vaporization method [521 other methods. 

Depending on the community, specific interests, and targeted applica- 
tions, nanotubes are regarded as either single molecules or quasi-onc-dimensional 
crystals with translational periodicity along the tube axis. As there are an 
infinite number of ways of rolling a sheet into a cylinder (two of them are 
represented in Fig. (|6.4p ) the large variety of possible helical geometries, 
defining the tube chirality, provides a family of CNTs with different diam- 
eters and microscopic structures. Some properties of these nanotubes, such 
as the elastic ones, can be explained within a macroscopic model of a ho- 
mogeneous cylinder. Others depend crucially on the atomic configuration. 
For instance, the electronic and transport properties, are certainly among 
the most significant physical properties of CNTs, and crucially depend on 
the diameter and chirality. This dependence on the atomic configuration is 
quite unique in solid-state physics. 

CNTs can be either semi-metallic or semi-conducting, with a band gap 
varying from zero to a few tenths of an eV, depending on their diameter and 
chirality. Further, the band gap of semi-conducting tubes, or the energy 
difference between the peaks in the electronic density of states, the so-called 
Van Hove singularities, can be shown to first order to be simply related to 
the tube diameter. Such remarkable results can be obtained from a variety 
of considerations, starting from the so-called Zone Folding Approach (ZFA), 
based on knowledge of the electronic properties of the graphene (a single 
sheet of graphite), to the direct study of nanotubes using semi-empirical 
Tight-Binding (TB) approaches. The comparison with more sophisticated 
Ab-initio calculations, and with available experimental results, permits to 
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Figure 6.5: The Zone Folding Approximation. On the left panel the section 
of graphene which describes the (3, 5) CNT within the ZFA. Periodic boundary 
conditions are imposed on the dotted lines. On the right panel the Brillouin Zone 
of graphene with the fc-points lines which respect the the rolling condition of a 
(5,5) CNT. From Ref. [2]. 

find the limits of these simple approaches. 
The Zone Folding Approach 

In the ZFA the wave-functions v]>„ ^.^ j,^ (x) = e*('^="^+''^^)M„_fc^^fe^ (x) of the 
graphene sheet are used to describe electrons in CNTs assuming that the cur- 
vature of the sheet gives negligible effects if the tube radius is large enough. 
Then the only difference the electrons feel, with respect to graphene, is the 
quantization of the angular momentum Iz = k^R (we are assuming here 
that the tube direction is the z-axis). Accordingly only the wave-functions 
gi(j^/fl<>+fc^z)y^^^^^^^^^ ^jj.^ are considered in the model. 

In Fig. (|6.5p . on the left (a) panel, a CNT is represented as a stripe in the 
plane of graphene. In the ZFA each CNT can be identified by two numbers, 
which represent the circumference vector in the basis of the direct lattice 
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vectors ai,a2 of graphcne. Boundary conditions are tlren imposed wliicli 
results in selecting specific fc-points in the reciprocal space. In Fig. (j6.5p . on 
the right top panel, the Brillouin Zone (BZ) of graphene is represented and 
the lines correspond to the k points which respects the boundary conditions 
of a (5, 5) CNT. The right bottom panel, instead, shows the energy surface 
of the TT — TT* bands of graphene, cut by the allowed /c-points lines which can 
be used to construct the CNT band structure. 

The predictions of the ZFA give a good description of the properties 
of CNTs but some corrections have to be considered both for SWCNTs, 
to include the curvature of the tubes, and for MWCNTS, to include the 
effect of the interaction among different tubes. For example the position 
of the so called Dirac points (see the next subsection for the definition of 
the Dirac points) has to be shifted in the Brillouine Zone due to curvature 
effects. Finally the ZFA results depend on the method used to compute the 
band structure of graphene. The energy surfaces in the Fig. ()6.5p have been 
calculated within the TB model, for example. 

The Tight Binding Model for graphene 

The graphene plane is an hexagonal lattice with two atoms per unit cell 
(A and B) and a basis defined by the vectors (01,02), as in Fig. 16.51 The 
condition a^ • hj = 2Sij allows one to obtain the reciprocal lattice vectors 
(61,52). Every carbon atom possesses four valence electrons (two 2s and 
two 2p electrons). When the atoms are placed onto the graphene hexagonal 
lattice the electronic wave functions from different atoms overlap. However, 
such an overlap between the pz orbitals and the s or and Py electrons 
is strictly zero by symmetry. Consequently, the Pz electrons, which form 
the TT bonds in graphene, can be treated independently of other valence 
electrons. Within this 7r-band approximation, the A atom (/ B atom) is 
uniquely defined by one orbital per atom site pz{r — va) [/ Pzir — tb)]- To 
derive the electronic spectrum of the total Hamiltonian, the corresponding 
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Wave vector 

Figure 6.6: Band structure of graphene evaluated with tiie tigiit binding model 
compared with the result of ab-initio calculations. From Ref. [2]. 



Schrodinger equation has to be solved, and by applying the Bloch theorem, 
the wave functions can be written as follows 

vl/k(r) = CA(r)p^,k(r) + CB(r)pf^k(r), (6.7) 

where 

pf,k(r) = -^^e^'^V.(r-rj-l) J ^ A,B . (6.8) 

k is the crystal momentum, TV is the number of unit cells in the graphene 
sheet, and 1 is the cell position index. 

The spectrum is derived by solving the Scrodinger equation which reduces 
to the diagonalization of a 2 x 2 matrix 



/ Haa-E Hab \ 
\ Hba Hbb ^ E J 



(6.9) 
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in the space defined by the |J) = P^|j.(r) wave-functions. Neglecting the 
overlap {A\B), restricting interactions to nearest neighbors only and setting 
Haa = Hb b = as energy reference the dispersion relation are then: 

2cos(ka2) + 2cos(k[a2 — ai])- (6.10) 

These are the tt and the tt* bands in the TB model. One of the two bands, 
which represent the valence and the conduction bands, is completely filled 
and the other completely empty. Moreover they intersect only in two points, 
known as Dirac points, in the BZ, as shown in Fig. (j6.5p . For this reason 
graphene is a semi-metal with a one-dimensional Fermi surface. For the 
same reason CNTs are predicted to be either metallic or semi-conducting 
according to whether or not the Dirac points belong to the set of fc-points 
allowed in the ZFA. 

The predictions of the TB model are partially confirmed by Ab-initio 
calculations even if some differences appear. For example the DFT band 
structure of graphene is not symmetric with respect to the chemical poten- 
tial, as the TB one (Fig. (j6.6|) ). In the rest of this chapter we will work in 
the ZFA using as starting point the TB band structure. Ab-initio corrections 
will be discussed later. 

Within the TB+ZFA scheme, it can be show that for example all (n, n) 
and (3ri, 0) CNTs are metallic, with n any integer, while all the remaining 
(n, 0) tubes are semi-conducting. CNTs of different kinds, (n, m) CNTs, are 
said to be chiral and can be either metallic or semi-conducting. The general 
rule is that a CNT is metallic if n — to is a multiple of 3 [5] . In Fig. (j6.7p 
the bands structures for the (5, 5), the (9, 0), the (10, 0) and the (8, 2) CNTs 
are shown. 
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6.3 Theoretical predictions and experimental 
results 

The Zone Folding Approach 

The state-of-the-art theoretical approach to the AB osciUations in CNTs is 
based on the ZFA. 

The ZFA is introduced by observing that the Hamiltonian of a CNT can 
be expressed in cyhndrical coordinates as 
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Figure 6.8: Aharonov-Bohm effect witliin the Zone Folding Approacfi. From 
Ref. [g. 



Then in the ZFA the foUowing map is apphed 



(j) -J> x/RcNT 

r 

z ^ z 



(6.12) 
(6.13) 
(6.14) 



where {x,y, z) are cartesian coordinates witli tlie graphene sheet oriented in 
the xz~plane. The ZFA Hamiltonian, when a magnetic field is present, is 
then 



H 



2m 



d 



dx ^oRcNT 



dz"^ 



+ V{x,y,z). (6.15) 



As the ZFA approximates the CNT with a planar graphene sheet in Eq. 
(|05l) the 2"^* term of Eq. (fCTTj) is set to zero and r w Rcnt- Apply- 
ing the resulting Hamiltonian to the block wave-function ^fc^^^^ (a;, y, z) = 
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(a) Metallic CNT (b) Semiconducting CNT 



Figure 6.9: Aharonov-Bohm gap oscillations in carbon nanotubes according to 
the zone folding approach. 



e 



"fc^.fc, (x) and defining 



K = kx- (6.16) 

we obtain, for the periodic part of the wave-function, an Hamiltonian iden- 
tical to the case without magnetic field a part from a shift of . 

This means that the eigen-functions of the Hamiltonian with B ^ can 
be written in term of the eigen-function of the Hamiltonian with B = 0: 

<..,fe.(^,y,^) = e'('===-+^-=^)<,^,,^(x) (6.17) 



e 



"^.I^,(x)e'(*/*°)(^/^^ (6.18) 



Consequently we can obtain the eigen-functions and the corresponding shift 
of the energies of the allowed fc-point grid as shown in Fig. (|6.8p . This means 
that, increasing the magnetic field, all the electronic properties of the CNTs 
will oscillate with period $o as, each times <I> = n$o the allowed fc-points will 
coincide with the ones at B = 0. hi particular the electronic gap is predicted 
to oscillate because periodically k'^ match a Dirac point. For 0) or (n, n) 
CNTs, two kinds of oscillations can exist according to whether the CNT is 
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(a) Experimental setup 



(b) Resistivity oscillations 



Figure 6.10: Description of tlie experimental setup and results of the experiment 
carried out in [66]. In tiie rigiit panel resistivity oscillations witii a period $o/2 
can be observed. 

metallic or semi-conducting as shown in Fig. (|6.9|) 
The experimental evidences 

The first experimental evidences of the AB effect in CNTs is reported in Fig. 
(|6.10p . In their work A. Bachtold et al. [SS] showed that the resistivity of 
a MWCNTs oscillates when an increasing magnetic field is applied. They 
measured oscillations with a period of h/2e in agreement with the prediction 
of the AAS theory. Assuming that the current carriers where localized, in 
the radial direction, on a single shell, they estimated their average radius 
from the relation nr^B = $o/2 = h/2e, valid after a complete oscillation. 
Knowing the external magnetic field they obtained an average radius for the 
current carriers of r = 8.6 ± 0.1 nm, in excellent agreement with Atomic 
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Force Microscopy measurements of the most external tube in the MWCNT: 
r = 8 ± 0.8 nm. That is, assuming that the carriers were mainly localized 
on the external tube, the AAS oscillations were confirmed. 

The results obtained by A. Bachtold et al. [66] however where not an 
experimental evidence of the gap oscillations existence in CNTs under the 
effect of an increasing magnetic field. Indeed the AAS effect is independent 
on the electronic gap and AAS oscillations can be measured in mesoscopic 
metallic rings too, contrary to pure AB gap oscillations which are a peculiar 
characteristic of CNTs. 

The first experiments able to test the predictions of the ZFA for the 
gap oscillations were performed in 2004. Coskum et al. [51] prepared an 
experimental setup (see Fig. (|6.11l o)) similar to the one of Ref. [66] but 
using MWCNTs shorter than the dephasing length to probe a qualitatively 
different phenomenon, namely, the electronic energy spectrum modulation 
by a coaxial magnetic field. In particular they measured the differential 
conductance at finite bias on a single-electron tunneling transistor formed 
by a MWCNT acting as coherent Coulomb island. In this way they were able 
to observe the interconversion of semi-conducting and metallic nanotubes. 
The results are shown in Fig. (|6.11l &'). where the differential conductance 
is plotted in units of e'^/h. Oscillations of the conductance with the correct 
AB period, $o = h/e are clearly shown in the figure. Moreover they were 
able to observe the splitting due to the interaction of the external field with 
the electron spin B • cr. 

Zaric et al. instead |70| measured the photoluminescencie spectrum of a 
SWCNT immersed in a magnetic field. They showed that the lowest energy 
excitonic peak (peak ^^1 in Fig. (jG.llI c)) moves to lower energy when the 
magnetic field increase. This result is in agreement with the prediction of 
the ZFA and in particular the gap closing in semi-conducting CNTs which, 
at low magnetic fluxes, is predicted to follow the rule Eg = Eo{l — 3$/$o)- 
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(a) Experimental setup (b) Conductance 




0.95 1.00 1.05 1.10 1.15 1.20 1.25 
Energy (eV) 



(c) Photoluminescencie spectrum of a SW-CNT 

Figure 6.11: Top panels: description of the experimental setup and results of the 
experiment carried out in In panel (b) the differential conductance oscillates 
with the imposed magnetic field in a standard circuit with an applied bias poten- 
tial at two different values of the gate potential. Low panel: photoluminescencie 
spectrum of a SW-CNT immersed in a static magnetic field from Ref. [7D]. The 
beginning of a gap oscillation is displayed by the behavior of peak #1 
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Chapter 7 

Numerical results 



7.1 Details of the implementation 

In order to describe magnetic field effects "ab~inito" we used tlie Yambo [HI] 
code. Yambo starts from a previous SCF computation of the ground state 
of the system at zero magnetic field, taking as input the KS (LDA) wave- 
functions Tpi and energies e^. Then it constructs a new Hamiltonian, H = 
HuFT + H„, where HnFT is the DFT Hamihonian with no external 
fields, and H„iagn = A"^^* • j The total Hamiltonian is constructed in the 
space of the KS wave-functions as 

H,^ = + V,f[5p] + (l^f '"^'"[p] - V^f '""[po]) + i?™'^^", (7.1) 

with po the density at B = 0, p the self-consistent density and Sp ~ p — po. 
At a first step p = po, then the Hamiltonian is diagonalized and a new set 
of wave-functions is obtained, together with a new p. A new Hamiltonian 
is then constructed from the new density and projected in the space of the 

^ Other considerations on the Hamiltonian in presence of a magnetic field can be found 
in Appendix IDI 



134 



Numerical results 



new set of wave-functions. The process is carried on until the convergence is 
reached. A convergence threshold for both the density end the eigen-energies 
is defined within the code. 

The implementation is done assuming that a small number of wave- 
functions in the KS basis-set are needed to have a good description of the 
system. Otherwise this implementation would be highly impractical as the 
code needs to store all the wave-functions of the basis set in memory and 
to compute all integrals numerically. Indeed, as long as the perturbation is 
small, the density of the system is not expected to change much and so a 
small basis-set is required. 

Among other reasons, our motivation to investigate ab-initio the effects 
of magnetic fields in CNTs is to test the correctness of the assumptions 
embodied in the ZFA. In the previous chapter, indeed, we have shown that 
the pure AB effect has to be measured with an experimental setup where 
a confined magnetic field is present and electrons move in a region where 
the magnetic field is strictly zero. Indeed the prediction of the AB effect 
had a strong impact on the scientific community because it introduced the 
possibility of an effect which is non local in the magnetic field. However the 
experimental setup in which CNTs are always studied is in sharp contrast 
with this assumption, since CNTs are fully immersed in a constant magnetic 
field. 

We have designed the code to simulate this two different geometries. 
First, we simulate the pure AB effect, where electrons travel in a space 
where B=0. Second, we use a uniform field. We will refer to the first imple- 
mentation as confined geometry, because the magnetic field is confined inside 
the CNT and null outside. This setup will be compared with the standard 
experimental setup where the magnetic field is uniformly distributed (we will 
refer to this case as extended geometry) . 



7.2 Gap oscillations 



135 



7.2 Gap oscillations 

We consider five CNTs: two metallic, the (5,5) and the (8,8), two semi- 
conductive, the (8, 0) and the (14, 0) and one Muhi-Walled, the (5, 5)@(10, 10) 
CNTs ground state have first been computed with the Abinit [Sll code us- 
ing a super-cell with dimensions (10 A.+2Rcnt) x (10 K+2Rcnt) x hcNT 
and an angle of 120° between the vectors defining the super-cell in the xy 
plane in order to maximazc the distance between CNTs in the periodic ar- 
ray of super-cells. Here Rcnt is the radius of the CNT while hcNT is the 
dimension of the supercell in the periodic direction, z. hcNT ~ 2.46 for 
metallic (n, n) CNTs where we used a reciprocal space grid of 20fc-points 
in the Brillouin zone, and hcNT = 4.26 for semi-conducting {n, 0) CNTs 
where we used a grid of lOfc-points. However finer grids and bigger cells in 
the xy directions have been tested when needed. Self Consistent calculations 
have then been performed with the Yambo code using a basis-set of about 
40 states per fc-point, above the last occupied state. 



Metallic single wall carbon nanotubes 

First we consider two metallic CNTs: a (5,5) tube with radius 3.39 A and a 
(8,8) tube with radius 5.41 A. In Fig. (|7.ip we compare the gap dependence 
on the applied magnetic flux in the two geometries with the result of the 
ZFA. In the case of the smaller (5,5) tube we immediately see a first impor- 
tant difference between the extended geometry and the confined geometry. 
The extended geometry, which represents the standard experimental setup, 
overestimates by ^7% the elementary flux $o which defines the periodicity 
of the gap oscillations. 

To explain the different gap dependence obtained in the two geometries 
we introduce, in a formal manner, the Hamiltonian which governs the AB 
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Figure 7.1: Gap oscillations in the metallic (5,5) and (8,8) CNTs. Two different ge- 
ometries are considered. In tlie extended geometry tfie magnetic field is applied uniformly 
in all the space (blue boxes). In the confined geometry, instead, the magnetic field is 
confined inside the CNT (red spheres). We compare the ab-initio calculations with the 
ZFA results (black line). We see that in the extended geometry, which represents the 
standard experimental setup, the Lorentz correction (see text) induces an overestimation 
of the elementary magnetic flux <I>o- 



effect in the specific case of a CNT: 



Eq. (|7.2|) describes tfie electronic dynamics under the action of a static mag- 
netic field, written in cylindrical coordinates centred in the center of the 
CNT. A (r) is the vector potential which, in the symmetric gauge, describes 
a static magnetic field along the z direction and V{r, ip, z) is the local DFT 
potential, which includes the ionic potential plus the Hartree and exchange- 
correlation terms. 

The only term of Eq. ()7.2p which reflects the different geometry [extended 
or confined) is A (r) . In the extended geometry 




(7.2) 



A 



Extended 



(7.3) 
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Figure 7.2: Two dimensional plot of the last occupied and the first unoccupied band 
wave— function at the F point for the (8, 0) CNT. The tt orbitals, deformed by the cur- 
vature of the CNT, have a larger amplitude in the outer part of the CNT surface. As a 
consequence we have a finite difference between the Hamiltonians corresponding to the 
extended geometry compared to the confined geometry. 



with i3o = |B|. In the confined geometry, instead, we have 



with $ = ttBqR'^j^j, and Rcnt the CNT radius. From Eqs. (fO)) and (fTi)) 
we see that A^^*^"''^'^ {Rcnt) ^ A^°''f''""^ (Rcnt), which imphes that, if 
the electrons would exactly move on the tube surface the extended geometry 
and the confined geometry would lead to the same gap oscillations. The 
different gap dependence observed in Fig. (|7.ip . is then due to the different 
correction induced in the total Hamiltonian by ^-Extended ^^^^^ ^Confmed^ 
If we plug the two different expressions for A (r) into Eq. (j7.2p we get two 



■^The potential in Eq. 1)7. 4p is not defined at r = 0. However this problem can be 
easily overcome in the implementation, sotting the magnetic potential A = l/2{^/rQ)r 
for r < ro with rg a tunable parameter. If ro is small enough then the vector potential is 
different from the AB potential only in a small region where the wave-function is almost 
zero and so the results are independent of ro. In our calculation we checked rg = 0.5Bohr 
is enough. 




1 $ 



(7.4) 



138 



Numerical results 



different Hamiltoiiians, ijcon/ined ^^^^^ ^extended ^ wliose difference is 

A U TTextended Tjconfined ekBp ( d \ ( ^CNt \ 

AH = H -H =-__^^2^---j (^1-— ^j, 

(7.5) 

with $0 = h/e. This term is zero when r ~ Rcnt, while near the tube 
surface behaves like ~ Bq (^l — ^'^"^ ^ ■ Now, as shown in Fig. (|7.2p . the 
electrons are localized near, but not exactly on the tube surface. Conse- 
quently 1 - {RcMr/r^) ^ and, (i/^^i^tended _ jjconf^ed^ ^ 

We will refer to the correction defined by Eq. (|7.5p as Lorentz Correc- 
tion (LC) as it introduces a magnetic term which depends on the electronic 
trajectory (through the term -r^eajtended-j^ rpj^^ appears in the extended 
geometry as an effective different radius of the electronic orbitals, as the cor- 
rection would be zero defining the flux with respect to the effective radius 
which satisfy the equation 1 — {R'^j^/r^) = 0. In the confined geometry this 
correction is zero by definition; physically this is related to the fact that the 
AB effect docs not depend on the specific electronic trajectory. 

The LC goes to zero in the limit of SW-CNTs with infinitely large radius. 
Nevertheless, even in large CNTs, impurities or defects can alter the elec- 
tronic trajectory creating deviations from a perfect circle of radius Rcnt- In 
all this cases we predict the LC to induce deviations from the AB oscillations 
when uniform B field is applied. 

From Fig. ()7.ip we see that the ZFA matches the ab-initio simulation 
of the pure AB effect corresponding to the confined geometry setup. This 
agreement is due to the fact that in the ZFA the LC is strictly zero as the 
electrons are assumed to move exactly on the graphite sheet, i.e. on the 
CNT surface. Consequently, in the ZFA the electronic gap is function of the 
flux only. 
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Figure 7.3: Gap oscillations in a semiconducting CNTs: (8,0) and (14,0). In 
contrast to the metallic case the curvature effects induce more evident differences 
betweem confined and the extended geometries. The same conventions of Fig. 17.11 
are used here. 



Semi— conducting single wall carbon nanotubes 

We now consider two semi-conducting CNTs: the (8,0) and the (14,0). The 
flux dependent electronic gap is shown in Fig. (|7.3p . Similarly to the metallic 
case, LC makes the extended geometry to oscillate with a period greater than 
h/e. In contrast to the metallic case, the gap vanishes at two values of $, 
which the ZFA predicts to be at $o/3 and 24>o/3, when the Dirac points 
become allowed k points Noticeably both points are renormalized in 
the ab-initio simulation by curvature effects. It is well known, indeed, that, 
compared to graphene, curvature effects shift the Dirac points |2| K (see 
also Fig. (|7.8p ) at a position \K\ < 27r/3a, with 27r/3a being the Dirac point 
position in graphene. Accordingly a lower magnetic field is needed to bring 
the Dirac point in coincidence with the set of the allowed fc-points, and 
semi-conducting CNTs become metallic at $ < $o/3. Being the oscillations 
symmetric, the second metalization point is reached at $ > 2$o/3. 

The deformation of the oscillations in the (n, 0) CNTs is smaller in big- 
ger tubes. However it goes to zero slowly, because both the shift and the 
magnetic period depend on the size of the tube. For this reason the effect is 
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still not negligible in the large (14, 0) tube, as shown in Fig. ()7.3p . 

From Fig. ()7.3p we see that the (8, 0) gap oscillations strongly deviate 
from the ZFA that does not reproduce, even qualitatively, the full ab-initio 
results. The reason for this large discrepancy traces back to the presence of 
a metallic-like band located near the Fermi surface. 

This band is shown in Fig. (|7.4p together with the tt/tt* bands closest to 
the Fermi level. When the B field is increased we see that, in contrast to the 
tt/tt* bands, the metallic-like band does not shift, but moves inside the gap 
produced by the tt/tt* states. Consequently by changing the flux intensity 
the gap is defined by transitions between the tt/tt* states or between the tt 
and the metallic-like band. This explains the anomalous dependence of Eg 
by in Fig. (fTlS]) . 

Multi wall carbon nanotubes 

Ahhough SW-CNTs arc routinely syntetized, MW-CNTs stiU constitute 
the majority of cases used in the experiments. In Fig. (|7.5p we consider 
the case of a (5, 5)@(10, 10) CNT with radn 3.39 and 6.78 A. In this case 
the confined geometry is implemented considering a fiux $ = ttBqR^-^^q ^^^^ . 
This flux is roughly the same experienced, in the extended geometry, by the 
electrons of the (10, 10) CNT. In the case of MWCNT the ZFA is used only 
by considering two different Hamiltonians, one for the (5, 5) and one for the 
(10, 10) and by neglecting the tube-tube interaction. Our ab-initio results, 
instead, reveal a quite different picture. 

The gap calculated in the extended geometry follows the ZFA prediction 
except in the very low field regime and near the first inversion point ($ w 
$o)- While the LC causes the shift of the inversion point, the metallic 
regime observed for $ < $o/10 is not described at all by the ZFA. It is, 
indeed, a consequence of the different chemical potential felt by the electrons 
moving on the (5, 5) and the (10, 10) surfaces. Here, moreover the LC is 
enhanced because the presence of the inner CNT tends to attract electrons 
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Figure 7.4: Metallic-like (in red) and tt / tt* (in blue) bands of the (8, 0) CNT. 
In presence of a magnetic field electrons spinning clockwise around the CNT have 
a different energy with respect to electrons spinning counter-clockwise. On the 
left panel $ = 0, while on the right $ > 0. We notice that the upper metallic- 
like band lays in the gap without acquiring a splitting due to the breaking of the 
time reversal symmetry, i.e. the energy of electrons spinning around the tube in 
opposite direction is different when a magnetic field is present. 



and accordingly to modify the effective radius defined by the LC. Indeed the 
flux renormalization is « 4% here, where i?(io,io) — 6.78 A, while in the SW 
configuration the renormalization is sa 7% and w 1% for the (5, 5) and (8, 8), 
where i?(5^5) = 3.39 and R(s,8) = 5-42 A. 

The ab-initio calculations show that there is a w 0.09eV shift between 
the two chemical potentials, as shown in Fig. (|7.6p . Here the band structure 
of the double walled CNT is plotted at $ = 0, on the left, and at small flux on 
the right. At zero magnetic fiux there are two pairs of crossing bands which 
can be identified as the tt — tt* bands of the (5, 5) and the (10, 10) CNTs 
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Figure 7.5: Gap oscillations in a Multi-Walled CNT: (5,5)@(10,10). In contrast 
to the SW-CNTs here the confined and the extended geometry display a qualitative 
different behavior. In the inset a complete period of oscillation according to the 
ZFA which describes the MW-CNT as two not interacting SW-CNTs. The same 
convention of Fig. 17.11 is used here. 



respectively. The two crossing points hovifever are not aligned in energy so 
that vifhen a small magnetic flux is present two small direct gap opens but 
the CNT remains metallic as long as the tip of the tt* band of the (5, 5) is 
lower in energy than the one of the tt band of the (10, 10) CNT. Even when 
a gap opens the system remains an indirect gap semiconductor for certain 
range of values of the applied magnetic flux. 

In the confined geometry the gap dependence is deeply modified. In the 
low $ regime a metallic region appear for the same reason outlined above. 
At difference with the extended geometry case, however, the gap increases 
with increasing $ with a much larger slope, dominated by the (10, 10) CNT. 
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Figure 7.6: Band structure of the metallic (5, 5)@(10, 10) CNT at zero magnetic 
flux (left panel) and at low magnetic flux, extended geometry (right panel). In the 
left panel in blue the bands, near the Fermi level, occupied from electrons on the 
inner shell, while in red the ones occupied from electrons on the outer shell. When 
the magnetic field is present it is not possible anymore to distinguish to which shell 
belongs each band. Due to the shift in energy between the point where the tt-tt* 
bands of the two CNTs cross, at low flux, though two small direct gaps open, the 
system remains metallic until an indirect gap opens (right panel). Only at higher 
fluxes the system becomes again direct gap semi-conductor. 

The slower slope of the extended geometry, instead, is dictated by the gap 
of the (5, 5) CNT. 

Another drastic difference with the extended geometry case is a second 
metallic phase near $ « This phase is due to the fact that, in the 
confined geometry, the electrons on the two CNTs feel the same magnetic 
flux. As a consequence the gaps of the two CNTs move coherently as <I> ap- 
proaches $0 and the same situation of $ « occurs. The coherent variation 
of the gaps of the (5, 5) and the (10, 10) tube result, for the gap oscillations. 
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in a period which is 1 /4 the period of the extended geometry. 



Numerical instabilities: a "gauge fixing" solution 

The results from the previous section have been obtained with an "improved 
implementation" which takes advantage of the gauge freedom (see Appendix 
[Ujl to describe a magnetic field. 

We observe that the AB effect can be described (as we did in the confined 
geometry) by a vector potential 

ft $ 

- --^ (7.6) 

which is a pure gauge field when $ = n$o with n any integer. A potential 
is pure gauge if a function A(r, -0,^) exist, which plugged in Eqs. (jC.f |) 
nullifies it everywher(|§. In this situation we need a function which satisfies 
the relation l/r dA/dcj) -A^{r) with defined by Eq. ((7^ . 

Indeed such function can be constructed for any value of the magnetic 
flux 

A(r>,z) = -5^</). (7.7) 

However in quantum mechanics a gauge transformation in the electro-magnetic 
potentials has always to be realized together with the corresponding trans- 
formation of the wave-function. The new wave-function obtained using Eqs. 
(|7.7p is well defined only if $ = n^o, while for any other value of $ it is a 
multi-valued wave-function. For this reason any effect induced by the AB 
effect is periodic with period $o- 

When a vector potential defined by Eq. ()7.6p with $ = 7i$o is applied, 
all the physical quantities of the system have to remain unchanged, while all 
wave-functions acquire a phase factor. However the KS basis-set is not a 



■^The scalar potential V must be zero. This condition is automatically satisfied in the 
static case, as a time independent gauge generating function can be used. 
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convenient choice to describe this change in the phase of the wave-functions. 

This can be understood if we try to construct the wave-functions at 
$ = $0 as a hncar combination of the KS wave-functions at $ = 0. The 
wave function of any single particle state must be expressed as 

If few elements of the basis-set are used to build i ^^"^ using a basis set 
which contains V'* 1^ and few other orthogonal wave-functions orthogonal, 
this leads to severe numerical instabilities that we first encountered when we 
tryed to compute the CNTs gap oscillations. 

The metallic CNTs were not returning metallic when $ = $o- Moreover 
for the semiconducting CNTs we found that using a small basis-set the gap 
oscillations, as well as other physical quantities, were not periodic in the 
flux. 

To solve this problem we observed that at $ = ?i<I>o it is possible to impose 
a phase-matching solution. That is, instead of projecting the Hamiltonian, 
at $ = n$o, on the starting KS basis-set we can project it in a basis- 
set modified adding the needed phase factor. In this new basis-set the 
Hamiltonian will be diagonal. Then for intermediate values of the magnetic 
flux we can select the most appropriate basis-set according to which integer 
values of the fiux we are closer. This numerical procedure avoids the self- 
consistent cycle to remain trapped in the initial gauge. 

In Fig. ([777)1 the gap oscillations of the (14, 0) CNT with and without the 
phase fixing implementation are shown. Without fixing the gauge the gap 
oscillations seem to converge towards the gauge fixing solution increasing the 
number of bands. However the convergence is very slow as very few wave- 
function from the basis-set can be used for computational limits. Indeed if 
in the plane-wave basis-set thousands of states are needed to reproduce the 
KS wave-functions, we expect that thousands of states would be needed to 
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Gaiige--fixing .solution: 140 bands 
Standard approach: 140 bands 
Standard approach: 170 bands 
Standard approach: 200 bands 



Gauge fixing solution: 90 bands 
Standard approach: 90 bands 
Standard approach: 1 10 bands 
standard approach: 130 bands 




Figure 7.7: The gap oscillations of the (14,0) CNT (the left panel) and of the 
(5,5) CNT (right panel). The results of the gauge-fixing implementation are com- 
pared with the results obtained with the standard implementation. Increasing the 
number of KS wave-function in the basis-set the two methods converge to the 
same result. 



correctly reproduce a phase factor. 

The same is true for metallic CNTs. For these the convergence problem 
in the standard implementation appears for field values corresponding to the 
metallic phase, as shown in Fig. ()7.7|) . 



7.3 The band structure 

The electronic gap of a CNT depends only on the behaviour of the last 
occupied and the first unoccupied band. In particular the gap depend of 
their behaviour near a specific /c-point, that is near the Dirac point (see 
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Figure 7.8: Band structure of the (5,5) CNT at zero magnetic field near tlie 
Fermi level. The Dirac point, according to the prediction of the ZF model, should 
be placed where the two red lines cross. Due to curvature effects, not included in 
the ZF model, the Dirac point, i.e. the point where the n and the vr* bands cross, 
is shifted in the DFT approach. 



Sec. [Oil. 

We have already seen how some deviation of the gap oscillations from 
the ZF predictions, can be understood in terms of the bands structure, as in 
the case of the (8,0) CNT or of the muhi-walled (5, 5)@(10, 10) CNT. 

In this section we will explore more in details the effects of the magnetic 
field on the bands structure CNTs. As a reference we will use the prediction 
of the ZF model. In Fig. (|7.8|) for example we see that the band structure 
of the (5,5) CNT computed with Yambo is similar to the prediction of the 
ZF model (see Fig. (6.7(a) ) but with the shift of the Dirac K point. This is 
a result known from the literature. 

When a magnetic field is present, then the band structure follow, in 
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Figure 7.9: Band structure of the (5,5) CNT near the fermi level at increasing 
values (black-pink-blue-violet-red-orange dots) of the magnetic field. Here the 
magnetic flux increase from <£> = (black dots) to $ = $o/2 (orange dots). 



first approximation, the predictions of the ZF model (see Fig. (|7.9|) ) and 
the deviation are the one we have already explored in the previous section 
describing the gap oscillations of different CNTs. 

In this section we will show (i) the role of the non-local part of the 
pseudo-potential for a correct description of the more bounds valence elec- 
trons and (ii) the effect of a not perfect alignment of the magnetic field with 
the CNT on the band structure. 



The non local part of the pseudo— potential 

When the term Hmagn = A • j is included in the Hamiltonian the non-local 
part of the pseudo potential must be changed accordingly (see App. |D]). In 
the present work we followed the strategy of Ref. [M] where the non local 
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Figure 7.10: KS eigenvalues at the Dirac point for the (8,8) CNT as a function 
of the magnetic flux. The results obtained using the standard pseudo-potential 
V^£'°{r,r') (black squares) and the corrected pseudo-potential y^^(r, r') (blue 
squares) are compared. Near the Fermi level (left panel) there are not appreciable 
differences, while for the more bound states (right panel) there is a difference of 
about 20%. 



pseudo-potential reads: 

y^Jr,r') = ^y#r°(r,r')e'/^/— (7.9) 

R 

In order to explore the effect of the expression appearing in Eq. 17.91 
we considered the behaviour of the KS eigenvalue at the Dirac point, not 
only for the last occupied and the first unoccupied band, but for all the 
eigenvalues. As for the electronic gap, all eigenvalues must be periodic with 
period $ = $0 and the oscillation symmetric with respect to $ = ^o/S- 

In Fig. (|7.10p we can see the different behaviour of the KS eigenvalues 
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induced by the correction to the non local part of the pseudo-potential, 
Eq. (|7.9p . In particular for the deeper valence states the change of the gauge 
at $ = $0 (see Sec. 17. 2[ subsection on the gauge-fixing solution) induces 
a discontinuity of the derivative when Vj^j^°{r,r') only is used. This is an 
indication that the KS eigenvalues are wrong in this case. The corrected 
pscudo-pocntial fixes the problem. 

The effect of a magnetic field not aligned to the tube axis 

Yambo [SI] is a plane-waves based code, devised to treat periodic systems. 
The applied uniform magnetic field is described with a vector potential of 
the form = Bpr (with r the radius in cylindrical coordinates) in the 
symmetric gauge. However such a vector potential is not periodic in the xy 
plane and numerically it is replaced by a saw-like dependence, which induces 
jumps at the super-cell borders in order to keep the potential periodic. 
CNTs are isolated systems in the xy plane and the electronic wave-function 
is almost zero on the borders, consequently the jumps do not affect the 
results. 

However, if the CNT is not perfectly aligned with the magnetic field, i.e. 
in the z direction a saw-like vector potential cannot be used. To overcome 
this problem it is possible to use an asymmetric gauge. Suppose for example 
that the CNT axis lays on the yz plane, then the vector potential 

Ay = Box (7.10) 

in the y-asymmetric gauge can be used. The generating function needed to 
switch from the symmetric to the y-asymmetric (/ x-asymmetric) gauge is 
A = ±1/2 Boxy . We have used this possibility in order to verify the effects 
of a possibly not perfect alignment between the magnetic field and the CNT 
in the experimental setup. First we checked that the different gauges gave 
the same results on a CNT in the z direction. Fig. (|7.1ip ). As expected 
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Figure 7.11: Convergence checks. Three different gauges describmg the same 
magnetic field with no appreciable differences. The test are performed on the 
(5, 5) CNT at low magnetic field in order to avoid wave-function phases problems. 



there are not appreciable differences between the three gauges. 

Then we considered a setup with a small misalignment between the mag- 
netic field and the CNT (6 = 15°), Fig. (f7T^ . In this configuration the 
small component of the magnetic field perpendicular to the axis of the tube 
has, in general, a negligible effect on the properties of the CNT. The band 
structure of the systems and the gap oscillations are dominated by the com- 
ponent of the magnetic field aligned with the CNT, that breaks the time- 
reversal symmetry and give the usual gap oscillation. However the perpen- 
dicular component is not always negligible, when two degenerate bands are 
present it breaks the accidental degeneracies in the band structure as shown 
in Fig. (Tmi) . 
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Figure 7.12: Low magnetic flux. A CNT perfectly aligned with the magnetic 
field in the left panel and the same CNT oriented with an angle of 15° at higher 
magnetic field on the right. The component of the magnetic field perpendicular 
to the CNT breaks the accidental degeneracies in the band structure. Here as 
example we show the effect on the last occupied bands. 

7.4 Persistent currents 

In the previous chapter we have introduced the concept of PCs induces by 
the AB effect. These currents have never been measured experimentally in 
SWCNTs because they are too small. However PC have been measured in 
mesoscopic rings [TSJ [771 HH] and have been predicted to be measurable in 
toroidal shaped CNTs jHS] within the TB model. 

We have compared the predictions of our ab-initio approach against TB 
results in order to estimate the reliability of the TB methods. Our scheme 
in fact has the advantage of including many-body effects and to describe the 
a states which are not in the TB model. 
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To compute PCs we've implemented two different scliemes: we evaluated 
the expectation value of the current-density operator j in order to construct 
the angular current and we also computed the angular current as the 
derivative of the total energy with respect to the applied magnetic flux. The 
second scheme includes many-body effect, that is the contribution due to 
the Hartree and the xc energies. The first scheme, instead, has been used to 
validate the method as the average on the KS wave-functions of the current 
operator can be compared with the contribution from the bands-energy, i.e. 
the sum of the KS eigenvalues. 

We will proceed as follow: (i) we will show how the total energy, and 
the different contributions to the energy, depend on applied the magnetic 
flux and how the limitations of our approach appear around 4> = ^>o/2. (ii) 
We will then show the results obtained for the total current using the two 
approaches. The preliminary results shown here arc for the metallic (8,8) 
CNT. 

The energy of the system 

The total energy of the system within DFT is expressed as [17] 

Etotip] = Stands - Eh[p] - {VM + EM (7.11) 

where Ebands = Eik ef '^(k)/i(k)u;(k) with ef '^(k) the KS eigenvalues, fi(k) 
the occupation factors and u)(k) the weights of the k-points in the BZ. 
The Hartree energy, Eh[p], is double-counted in the Elands term and is 
subtracted while the xc energy is incorrectly described by the V2;c[p] term 
which is subtracted. 

The different components of Etot [p] arc shown in Fig. I7.13[ that shows 
how the total energy (and also all other components of the energy) goes 
quadratically with $ and present a derivative discontinuity at $ = 'I'o/2. 
The quadratic behaviour can be understood, in first approximation, from 
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Figure 7.13: Energy components of the (8,8) CNT as a function of tlie applied 
flux in the confined geometry. From the top to the bottom: the total energy (black 
line), the Hartree energy (red line), the xc energy (green line), the sum of the KS 
eigenvalues (blue line), the expectation value of the xc potential (orange line). 



the structure of the Hamiltonian 

H,nagn - A • j = B • L + qpA^ . (7.12) 

At low magnetic fields the contribution from the term B • L is almost 
zero because all the contribution ~\~Lz and —Lz cancel almost cxactlj0 while 
the quadratic term grows. This behaviour, which is exact at low magnetic 
flux within the DFT scheme [l^, should be modified at higher magnetic 
fields. However in our approach the high fluxes correction to the quadratic 
behaviour is not correctly described, as witnessed by the discontinuity at 



*The wave— functions are not exact eigenstates of Lz because the rotation invariance is 
broken by the presence of the carbon atoms. 
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Figure 7.14: Current components of the (8,8) CNT as a function of tiie applied 
flux in the confined geometry. The expectation value of the j operator on the 
KS wavefunctions is compared against j — dEtands/d^ (left panel). The different 
components of the current are then computing deriving the different contributions 
to the total energy (right panel), cfr. Fig. 17.131 

<& = <J>o/2. As a consequence the gauge-fixing solution induces a unphysical 
jump in the current in order to restore the correct periodicity. 

The current 

The quadratic behaviour of the energy is reflected in the linear behaviour of 
the persistent current as a function of the magnetic flux (see Fig. 17.141 and 
the main components of the total current are the bands contribution and 
the Hartree contribution. First we can see that the two methods gives the 
same results for the bands term of the currenljfl. 

The result obtained both for the bands current and for the total current 
are much bigger than the prediction within the TB model |85| I76j . In par- 
ticular the self-induced flux, which within the TB scheme is around 10"'^$, 
is here close to 0.154>. This strong discrepancy can be due to many factors: 
in the TB scheme the many-body effect are not included and only the tt 
states are used to describe the total current. On the contrary in our scheme 

^Here it is plotted the total current in the unit-cell, that is the current for a length of 
2.46 Aof the infinite CNT. 
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the many-body effect are, at least in part, included and the sigma states 
are fully described, though we do not have an accurate description of the 
current in the region near $o/2. We are presently working on an improved 
implementation able to describe the current in the whole range of applied 
fluxes. 

7.5 Conclusions 

In this part of the thesis we have studied the effects of a magnetic field on the 
electronic properties of Carbon Nano-Tubes (CNTs) aligned with the field. 
In particular we focused our attention on the oscillations in the electronic 
gap induced by the Aharonov-Bohm (AB) effect. By using parameter-free 
approach we made a direct comparison between the pure AB interpretation 
(confined geometry) and the common experimental setup (extended geome- 
try). 

In the extended geometry we confirmed the general behaviour of the 
pure AB interpretation for ideally perfect Single Wall (SW) CNTs. However 
we revealed the existence of corrections due to a trajectory dependent term 
(the Lorentz Correction (LC)). This effect decreases increasing the dimension 
of the CNT in ideally defect-free CNTs, but it is likely to be enhanced if 
vacancies or impurities, which can alter the electronic trajectory, are present. 

In the literature the pure AB interpretation is usually used to describe 
magnetic field effects within the Zone Folding Approach (ZFA). It is known, 
however, that the ZFA does not take into account curvature effects that can 
modify the band structure of CNTs. We revealed these effects as corrections 
to the behaviour of the gap oscillations, both in the extended and in the 
confined geometry. The shift in the position of the metallization points in 
semi-conducting CNTs or the completely different shape of the gap oscilla- 
tions in the (8, 0) CNT arc examples of curvature effects. 

We also discussed how the interaction among different shells in a Multi 
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Wall (MW) CNT can alter the gap oscillations obtained in the SW con- 
figuration. Indeed in the MW configuration there is a drastic difference 
between the extended and the confined geometry. The AB interpretation in 
the extended geometry can be recovered only at the price of describing the 
MWCNT as two not interacting SW-CNTs. This is the standard approach 
used in the literature and again our approach has revealed important cor- 
rections like the existence of a metallic phase at low magnetic flux followed 
by an indirect gap phase. 

In the last two sections we have shown some preliminary results regarding 
the case of a not perfect aligned magnetic field and the existence of persistent 
currents in CNTs. Further work is planned in this direction in the near 
future. 

In conclusion we have described, for the first time at our knowledge, 
the AB effect in CNTs with a first principles approach. We showed how 
curvature effects modify the properties of CNTs under the effect of a mag- 
netic field. Our results are grounded to well-known facts in the ab-initio 
community. For example the shift of the metalization points we observed 
in semi-conducting CNTs can be related to the shift of the Dirac points in 
the band structure of small metallic CNTs. In the same way the peculiar 
behavior of gap oscillations we observed in the (8, 0) CNT is related to the 
existence of a metallic-like valence band, not predicted by the ZFA. 

On the other hand the existence of LCs, which can effect the electronic 
properties of CNTs is a completely new effect, which to our knowledge 
has never been addressed in the literature. Moreover in the case of the 
(5, 5)@(10, 10) CNT we pointed out how the pure AB interpretation of mag- 
netic field effects in MWCNTs is not free of ambiguities. Indeed the pure 
AB interpretation gives different predictions with respect to what is mea- 
sured experimentally in the extended geometry. The pure AB effect can be 
recovered, as a first order approximation, only under the assumption that 
electrons can be distinguished according to weather they orbit on the inner 
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or on the outer shell of the MWCNT. 
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Conclusions 



When I started to work on the subjects of the present thesis, more than 
three years ago, my idea was to carry on a project which would have been 
concluded at the end of the PhD. Now I realize that a research project, often, 
opens more question than the number of answers and solutions it founds. 
The work presented in this thesis is the result of a long learning process that 
is far from being closed. 

The initial plan of the thesis was to developed a new method able to 
describe double excitations, within a many body approach. This purpose 
required a deep study of the different approaches available in the literature: 
from the Bethe-Salpeter Equation and the time dependent density functional 
theory for extended systems, to the Hartree-Fock and Post Hartree-Fock ap- 
proaches to isolated systems. From the Configuration Interaction expansion 
and other wave-functions based approaches adopted in quantum chemistry 
to the random phase approximation and second random phase approxima- 
tion used to describe nuclear systems. 

Our resulting idea has been to merge some aspects of techniques used for 
extended systems with the key ingredients used in quantum chemistry and 
in nuclear physics. Indeed the two relevant aspects we focused on are the 
idea of screening and the mathematical properties imposed by exchange ef- 
fects. Following both physical intuition and mathematical rigor we proposed 
the resulting scheme as a possible choice to describe double excitations in 
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correlated materials. 

An interesting result is that the inclusion of double (and even higher 
order) excited configurations naturally emerge as a frequency dependent 
Bethe-Salpeter equation kernel in the space of single particle transition. 
This result is similar to the one obtained in other works [TT] were the fre- 
quency dependency, i.e. temporal non locality, appears from the contraction 
of an higher to a lower dimensional space. This projection implies that the 
system cannot be regarded anymore as closed, and so its evolution depends 
from its past history. 

The time-dependent density-functional theory and the Bethe-Salpeter 
equation kernels are frequency dependent because these schemes, instead of 
the many body wave-function, are based on the time dependent density and 
on the two particles Green's function respectively. Here we rediscovered the 
frequency dependence of the kernel as a consequence of the projection of the 
excited states operator in the space of single particles transitions. "There is 
a pleasure in recognizing old things from a new point of view. Also there are 
problems for which the new point of view offers a distinct advantage" |86| . 

Only systematic tests on realistic materials will reveal if the proposed 
approximation will work on realistic systems. This is an open question which 
will likely need much more time then the one available for a PhD thesis to find 
a definitive answer. When we arrived to formulate our final expression for 
the kernel of the Bethe-Salpeter equation we felt we had reached a satisfying 
point of our investigation. 

In the second part of the present thesis we tackled the description of 
magnetic field effects in carbon nanotubes aligned in the field direction, 
within a first principles approach. 

The state-of-the-art on the subject describes these effects in terms of the 
Aharonov-Bohm effect within the zone-folding Approach. Many of the prop- 
erties are obtained starting from a tight binding calculation on a graphene 
sheet, which can be performed analytically. We have shown how the first 
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principles description confirms the qualitative results of model. However 
corrections need to be considered. The first-principles approach includes, in 
a consistent manner, many effects beyond the zone-folding approach, such 
as the curvature of the graphene sheet in the carbon nanotubes or the cor- 
rection (Lorcntz correction) to the pure Aharonov-Bohm description of the 
magnetic field. 

In this part of the thesis the strategy adopted is, for certain aspects, 
opposite to the one adopted in the first part. Instead of mathematical rigor 
(the approximations used are, sometimes, not fully justified from a theoret- 
ical point of view) the key aspect here have been to focus on the physical 
behavior of the system. The description of a physical system often involves 
many different aspects and it is crucial to find out which of these are the 
most important. 

Moreover mathematical inspection has been used a posteriori to check 
and improve the approximation involved. This is the example of the "gauge 
fixing" solution which we have developed to overcome numerical instabilities, 
observing that the Aharonov-Bohm effect has to be a pure gauge effect for 
some values of the vector potential. 

The resulting approach enabled us to give a much accurate description 
of the system compared to what can be obtained using models such as the 
the zone-folding approach. 

This second part of the work has been an occasion to study a new subject 
and learn a different approach to a research project. The concept of double 
excitations had a clear mathematical definition though at the beginning we 
had no clue on how to incorporate them in a many body approach. 

On the other hand we had a clear idea of the experimental setup used 
to simulate carbon nanotubes immersed in a magnetic field, but it took 
some time to understand how the concept of the Aharonov-Bohm effect can 
be used to describe such setup. Indeed according to Aharonov and Bohm 
the Aharonov-Bohm effect arise because "there exist effects of potentials 
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on charged particles, even in the region where all the fields (and therefore 
the forces on the particles) vanish" [53]. However the experimental setup 
considered for carbon nanotubes is in sharp contrast with the situation sug- 
gested by Aharonov and Bohm themselves, since the carbon nanotubes are 
fully immersed in the magnetic field. Only understanding this difference we 
have been able to address specific questions and recognize, for example, the 
renormalization of the gap oscillations due to the Lorentz correction effect. 

Indeed when, at first, we decided to work on the description of magnetic 
field effects in carbon nanotubes our idea was to tackle the description of 
either magneto-optical spectra, and the "dark to bright transition" of the 
lowest energy exciton [57] due to the Aharonov-Bohm effect, or the resistiv- 
ity oscillations, related to the combined effect of Aharonov-Bohm and weak 
localization [551 E^. However, the description of the Aharonov-Bohm effect 
within a full ab-initio approach turned out to be an almost unexplored field, 
and we were forced to first address the many open questions which in the 
end became the core of the present work. 

As I stated at the beginning of this conclusions, many questions have 
been opened. I hope that the present thesis will be used as a starting point 
by someone, maybe myself, to look for new answers. Many possible paths 
have been highlighted, which could possibly inspire new projects. Started 
from the results presented in part H of the thesis a project apply the ap- 
proach proposed to describe double excitations on realistic materials. From 
part HI, on the other hand, it could be interesting to look for experimental 
configurations where the Lorentz correction is more pronounced than in ide- 
ally perfect carbon nanotubes. Similarly a new project could improve the 
ab-initio description of persistent currents, in order to say the last word on 
their real intensity and on the role played by many-body effects. 
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Appendix A 



Connection to the 
experiments: extended 
systems 

A.l The Dielectric constant 

In Ch. [T] we draw the connection between the microscopic quantities and 
the observable measured in an absorption experiment, that is the dielectric 
constant t{Lo) for extended systems and the polarizability a{uj) for isolated 
systems. Then in Ch.|3]we specialized the description to isolated systems and 
we wrote the expression of the polarizability in the space of single-particle 
wave-functions starting from the solution of Eq. ()4.8p . In sec. 15.41 we used 
this result to describe the spectra of two molecules {Gs,H2 and CiHo) and 
to test different aspects of the kernel proposed in Ch. [SJ 

However the kernel we propose, including the concept of screening, could 
in principle be used to describe extended systems too. For this reason we 
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introduce here the direct expression for the dielectric function obtained from 
the solution Eq. (|4.8p . It is convenient for this purpose to work in the space 
of the Block wave-functions: 

^„k(x) - «„k(x)e''^ - ^ c„(k + G)e('^+°)''. (A.l) 

G 

With this choice the structure k + G is introduced in the reciprocal space, 
and the four point response function can be expressed as 

^Gi,G2,G3,G4(ki,k2,k3,k4|a;) = ^ c„(ki + Gi)c*„(k2 + G2)x 

7im,st 

L„™,,t(ki, ka, kg, k4|cj)cs(k3 + G3)c*(k4 + G4). (A.2) 

Thanks to translation invariance we then reduce the dependence to three 
momenta k = (ka — ki), k' = (k4 — ks) and q = 1/2 [(ki + ka) — 
(k4 + ka)]. The contraction in real space to obtain the response func- 
tion x(l,2) = L(l,l;2,2) is equivalent to set (ki,Gi) — (k2,G2) and 
(ka, G3) = (k4, G4). So we obtain the function 

XGi,G2(q,'^) = LguGuG2,G2{(1,^)- (A. 3) 

We recall that the macroscopic dielectric function (sec Ch. [T|) is obtained 
by averaging over the unit-cell, that is, in this notation, taking the the 
component (Gi = G2 = 0) of the microscopic dielectric function. By using 
the expression 



Ci(ki + Gi) = 




(A.4) 



and rewriting the interaction in the space (k + G), wg,G' (q) = Sg,G' 47r/|k + Gp 



A. 2 Electron energy loss spectroscopy and absorption 
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we obtain 
eo,o(q''^) = 1 + ''«o(q)xo,o(q,w) 

(A.5) 

The absorption spectrum is described by the q —>■ limit of Eq. (|A.5|) . 
This is because the photon, at the range of energy of few eV, has negligible 
momentum pf = E/c ~ 10"^'' [kg m/s], if compared to the momentum of 
the electron at the same energy pe = y^2mei?) — 5 • 10"^^ [kg m/s], that 
is pf/pe^2-l0-^. 



A. 2 Electron energy loss spectroscopy and ab- 
sorption 

Once we computed the macroscopic dielectric function in extended systems 
electron energy loss spectroscopy (EELS) experiments can be described. In 
contrast to absorption experiments where the quantity Img[e\ is measured, 
EELS experiments measure the quantity Img[e^^] [7]- Here we show that 
this difference is related to how the long range term of the Coulomb in- 
teraction w enters in the equation to construct the macroscopic dielectric 
function. Indeed, using the microscopic relations e^^ ~ 1 + wx, Eq. (|2.33p . 
and e^/ = (e^i), Eq. (fTTOl) . we will derive for EELS 

<^Ai (q, = 1 + Wo (q)xoo (q, w) , (A.6) 
as opposite to the equation 



eA/(q,w) = 1 - wo(q)xoo(q,'^), (A.7) 
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which holds for absorption experiments. 

Here x is the response function obtained from the Dyson equation x = 
X° + X°(w + fxc)x where, with respect to Eq. p.31|) . the G = term of the 
interaction is not included in the Hartree part of the kernel, i.e. u'o(q) = 0. 

The derivation is a mathematical exercise, but the result has an intuitive 
physical interpretation. While in absorption experiments the photons probe 
the system locally, so they do not feel the effect of the long range contribution 
of the interaction, in EELS experiments the electrons travel through the 
medium for long distances and the system reacts in a different manner. This 
difference disappear when isolated system are considered, as long as the 
dimension of the system is smaller than the wave-length of the photons (i.e. 
around 10^ — 10^?ito). 

From Eq. (|A.6|) it is possible to understand why, sometimes, independent 
particle spectra are referred to as RPA spectra in the literature. At the RPA 
level fxc = 0, the term w of the kernel can then be divided in two parts: the 
long range contribution zi;o(q), which account for the long range part of the 
electron-electron interaction, and the other terms wa{ci) with G 7^ 0, the so 
called Local Fields (LFs), whose average on the unit cell is zero. If the LFs 
effects, which in extended system can be negligible, are not considered then 
X ~ at the RPA level. For this reason the IP approximation is sometimes 
referred to as RPA without LFs effects. 



Mathematical Derivation 

In order to obtain Eq. (jA.7[) we observe that the microscopic dielectric func- 
tion can be obtained inverting Eq. (|2.33p : 

2) = ^ = S{1, 2) - wil, l')n*(l', 2). (A.8) 



A. 2 Electron energy loss spectroscopy and absorption 
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We use here x, the retarded version of the T-ordered 11*, to write in recip- 
rocal space 

EG.G' = 5g,g> - WgXG,G'- (A.9) 
We use now the general rule for a matrix of dimension N x N 

(TO Ci \ , . . 

c. a)- 

where to = Mi^i, A is a matrix of dimension (N — 1) x {N — 1) and Ci and 
C2 are matrix of dimension 1 x {N — 1) and {N — 1) x 1 respectively, 



to express e^^ ~ l/^o starting from Eq. (|A.9|) : 

e^"^ = 1 - woXo.o - X! ^0X0, G £g]g' wg'Xg',0- (A.12) 

G,G'5^0 

Finally we recognize that 

Xo,o + J2 

Xo.G ^g]g' ■'^g'Xg',0 (A. 13) 

G.G't^O 

is the zero component of a modified Dyson-like equation where the G = 
component of the kernel is zero. Using the relation — (1 — wx)^^ in fact 
we can define 

X = X + XWX, (A.14) 
and so express the macroscopic dielectric function as 

e^^ = 1 - «;oXo,o. (A.15) 
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Appendix B 



On the quasiparticle 
concept 

To clarify the quasiparticle (QP) concept here we consider the poles of the 
one-particle Green's function (GF), that represent the description of photo- 
emission excitations. 

The QP concept introduced by Landau derives from the physical intuition 
that some excited states of an interacting many-electron system reassemble a 
one particle resonance in experiments carried out in accelerators. The main 
difference is that real particles do not interact with the background (the 
vacuum) while in the interacting systems QP excitations have a (long) finite 
life-time due to the interaction with the many-electrons sea. This means 
that QPs are not exact eigenstates of the Hamiltonian, but are "quasi- 
eigenstates" . 

In the GF formalism the finite lifetime is described by the imaginary part 
of the self-energy which, evaluated at a QP pole, is not real. This concept 
is strongly related to the existence of a continuum of poles in an infinite 
system, that is a branch cut in the complex plane. 
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First we give the mathematical construction which connects the Lehmann 
representation to the QP representation of a GF. Most of the concepts are 
in Ref. gS]. 

Consider the function (here x is in the complex plane) 

fA^) = (B-1) 

X — ei{k) + irj 

which has a pole on the real axis at a; = ei(fc) (the small if] off-set is there, 
as usual, for mathematical convenience). Now let's take the related function 

which has a series of poles on the real axis. Suppose that, in some limit 
(the thermodynamic limit in a physical system), the variable k becomes a 
continuum variable and the sum over k becomes an integraj^. Supposing 
for simplicity ei{k) = k 5i,i and that ai{k) ~ (l/7r)(l/[(fc — + the 
integrand has a branch cut in the lower complex plane for k ~ x — irj and 
a simple complex pole at k = Ei + iE2- We can perform a contour integral 
in the upper plane (avoiding the branch cut) and using the residue theorem 
we obtain 

TT U! — [El + lE2) 

where I let the 77 to go to zero as we do not need it any more here. 

The poles of g{x) can be either the branch cut, i.e. the poles of the 
Lehmann representation (Eq. (|B.2|) ). or the complex pole, i.e. the QP pole 
(Eq. (|B.3p ). This connects the Lehmann representation with the QP picture. 

We have a connection with the definition I gave in my thesis with Ei = 
Ci + i?e[Sii(i?i + iE2)] and E2 ~ Img[Y.ii{Ei + iE2)]- This last connection 



^In calculations on realistic system wc have a discrete grid of k points but we think 
this as an approximation to the whole Brillouin-zone 
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can be obtained writing the analytic continuation of 

w - (£i(k) + Eii[k,Lu)) 

The GF has a QP pole when Img{Y.ii{k,uj)) is small (wc need the QP to be 
almost an eigen-state of the Hamiltonian according to the request that we 
have a quasi-eigen-state). 

In the analytic continuation we find the complex poles and not the branch 
cuts on the real axis because the GF is defined on the Riemann surface where 
the branch cuts are the connections of two sheets. 

When we do not have a continuum of states, then all this procedure makes 
no sense and we remain with the simple poles on the real axis. Physically 
this is related to the fact that there is not a macroscopic number of states 
which can interact and make the QP poles to have a finite life-time. 

The definition of QP as "dressed" one particle excitations, which never 
consider satellites as QP, is not exactly the same as the idea of a quasi- 
eigenstate of the Hamiltonian. If a satellite is a clear and well definite peak 
in the spectrum, then it can be considered a QP. On the other hand if the 
satellite appear, for example, as a shoulder of a QP peak in the spectrum it 
cannot be considered a QP. 

However when, in practice, we do GW calculations we assume that the 
QP wave-functions are well described by the KS wave-functions. This means 
that we do not have the QP wave-function which describes the satellite. 
Maybe in this sense just "dressed" one particle excitations are QPs. 

The procedure I've described in this section is usually not carried on 
in the BSE scheme where instead the Lehmann representation is usec|§. It 
could be interesting to keep in mind this possibility, especially once one goes 
beyond the static approximation for the kernel of the BSE and so can think 
to an analytic continuation and to a kernel which can be non-real. 

^Only in some works the imaginary-part of the self-energy is used to evaluate the 
life-time of neutral excitations 
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Appendix C 

Gauge transformations 



In quantum mechanics a transformation is a gauge transformation if it leaves 
unchanged ah physical predictable quantities. The most well known kind of 
gauge transformation arc those related to the electromagnetic potentials (or 
more in general to gauge potentials). The Schrodinger equation, for example, 
is invariant under the transformations 



A ^ A + VA, (C.l) 
* ^ ^e^"^ . (C.3) 



V and A are the scalar and the vector electro-magnetic potential respec- 
tively, \1/ is the wave-function of the system, and A is a generic function. 

In Ch. [71 for example, we used this freedom to fix numerical instabilities 
in the solution of the self-consistent problem. The concept of AB effect 
itself is deeply related to the concept of gauge transformation. Indeed the 
existence of the AB effect is due to the fact that, in quantum mechanics, 
the vector potential outside a solenoid can be nullified only at integer values 
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of the magnetic flux trapped by the solenoid, and accordingly is not a pure 
gauge potential for non-integer values of the flux. 

The many-body Hamiltonian introduced in Ch. [T] is written approxi- 
mating the Hint term in the non-relativistic limit, i.e. neglecting terms 
proportional to 1/c (with c the speed of the light). The retained term, w, is 
the interaction in the "Coulomb gauge" or transverse gauge obtained setting 
V • A = and it is instantaneous. 

Other possible gauge transformations exist. The response function intro- 
duced in Ch.[T]in the dipole approximation is proportional to the expectation 
value (x) of the position operator. However the relation [H, x] = p between 
the position operator to the momentum operator can be used to write the 
response function in terms of either the position (lenght gauge) or the mo- 
mentum (velocity gauge). 



Appendix D 

DFT and magnetic fields 



In Ch. [7]we have shown ab-initio results of the effect of a static magnetic field 
on CNTs. The implementation in the code has been done adding the term 
Hmagn = A'^''* • j. While this term enters in the many-body Hamiltonian, 
it is not the only term which enters be used in the DFT Hamiltonian. Our 
implementation is then approximated. 



DFT vs CDFT 

When the coupling term A • j is considered the HK theorem does not hold 
anymore and CDFT have to be used, as we stated in Ch. [3l where an ex- 
tension of the HK theorem can be proved [HI [15] . Within (non relativistic) 
CDFT the total energy is a functional of the density and the paramagnetic 
current E[p,]^P'']. DFT can be seen as an approximation to CDFT with 
£;[pj(p) = 0]; accordingly Kc[p,.7^p^] ~ V^/^ and A^c ~ 0. This approxi- 
mation is used in the ZF model and in the present thesis too. 
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Diamagnetic screening 

Beyond the xc vector potential, the magnetic current generates a diamag- 
netic screening potential, that is the equivalent of the Hartree term which 
originates from the density: 

47r 7 k - I" I 

The diamagnetic screening term is usually negligible as the "diamagnetic 
field" is usually much smaller then the external applied field. In CNTs too 
is small, however, due to its cylindrical geometry it could be measurable. 
Indeed the "diamagnetic field" is, in first approximation, proportional to the 
radius of the electronic orbits. In simple geometries electrons spin around 
the atoms at a radius close to 1 Bohr. In CNTs and in mcsoscopic rings 
however the radius is much bigger. Indeed in mesoscopic rings the field can 
be strong enough to have stationary states with self-sustaining current^. 

In the present work we neglected this term as a first approximation, 
though the description of persistent currents in the ground state could be 
used to estimate the value of it. 

Non-local term of the pseudo— potential 

When a vector potential is considered, the non local part of the pseudo- 
potential has to be modified. This can be understood observing that a 
non local operator y„;(r,r') can be expressed as ]4;(r, p) [HHIIHI]- Then 
when a vector potential is present the substitution p ^ (p — eA) has to be 
performed in order to ensure that (y„;(r,p)) is gauge invariant and so that 
the total energy of the system is gauge invariant. 



^When the radius of the cylinder becomes very large and the classical limit is ap- 
proached the diamagnetic current must be balanced by the paramagnetic current and the 
total effect is not always diamagnetic. 
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We explored the effect of the non-local term of the pseudo-potential in 
Sec. [731 

The non -local term in the pseudo-potential enters even in the definition 
of the current operator. Indeed the current is proportional to the velocity 
and so to the commutator of the Hamiltonian with the position operator 



In the present work, when we compute the PCs we considered the effect due 
to the presence of the non-local term on the pseudo-potential. However we 
verified that this term is negligible with respect to the paramagnetic and the 
diamagnetic term of the total current. 



J = 



q ip-eA+ [Vni,v]). 



(D.2) 
(D.3) 
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